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Abstract 

Holloman  Air  Force  Base  houses  the  10  mile  long  Holloman  High  Speed  Test 
Track  in  Alamogordo,  New  Mexico  and  can  run  hypersonic  experiments  at  speeds  up  to 
10,000m/s.  Test  objects  are  loaded  onto  sleds,  which  are  connected  to  the  track  by 
“slippers”  that  slide  along  the  rails.  The  payload  sled  is  propelled  down  the  track  by  a 
series  of  rocket  sleds.  The  ability  to  predict  the  wear  that  will  occur  on  the  slippers  during 
these  experiments  is  important  for  slipper  design  and  preventing  catastrophic  failure. 
However,  high  speeds  complicate  wear  prediction  since  there  are  additional  contributing 
factors,  including  inconsistent  surface  contact,  fluctuating  thermal  and  friction  properties, 
and  additional  speed-induced  melt  wear.  The  goal  of  this  project  is  to  utilize  a  hydrocode 
program  -  CTH  -  to  predict  high  speed  wear  results  using  a  damage  criterion  that  allows 
damage  in  cells  to  occur  and  have  the  cell  removed.  Results  will  be  compared  to  the  wear 
seen  in  other  plane  strain  models  which  did  not  allow  the  damaged  regions  to  be 
removed.  Comparisons  are  made  to  a  3D  ABAQUS  simulation,  and  experimental  data 
derived  from  a  2008  test  run  at  the  Holloman  High  Speed  Test  Track. 
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PREDICTING  THE  WEAR  OF  HIGH  SPEED  ROCKET  SLEDS 


I.  Introduction 


1.1  Problem  Background 

The  Holloman  High  Speed  Test  Track  (HHSTT)  at  Holloman  Air  Force  Base 
(AFB)  in  New  Mexico  runs  hypervelocity  experiments  up  to  Mach  10  over  its  10  mile 
track.  These  experiments  provide  valuable  information  about  payload  performance 
inflight  environments.  However,  these  speeds  introduce  dynamic  wear  challenges  that 
have  not  been  thoroughly  investigated.  Figure  1  shows  an  example  experiment  on  the 
track. 


Figure  1:  Example  Experiment  Set-Up  at  the  Holloman  High  Speed  Test  Track  [1] 

The  setup  is  composed  of  multiple  sleds  mounted  on  a  steel  track.  The  front  sled 

has  the  experimental  payload,  and  the  trailing  sleds  have  rockets  loaded  on  them  that 

propel  the  fore  body  sled  down  the  track.  A  better  view  of  the  stages  is  seen  in  Figure  2. 
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Figure  2:  Stages  of  the  2008  HHSTT  Experiment  [2] 

The  connecting  piece  between  the  sled  and  the  track  is  the  slipper  and  is  shown 

below  in  Figure  3. 


Figure  3:  VascoMax300  Slipper  attached  to  the  1080  steel  Rail  [3] 

As  the  slipper  slides  along  the  rail  during  the  experiment,  the  bottom  surface 

progressively  wears  down  as  material  is  removed.  A  microscopic  view  of  the  slipper 

surface  provides  some  more  insight  into  this  phenomenon  in  Figure  4. 
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(a) 


(b) 


Figure  4:  Scanning  Electron  Microscope  photos  of  (a)  an  Unworn  Slipper  Sample  and  (b)  a  worn 

slipper  sample  [4] 

On  the  left  is  an  unworn  slipper  sample,  and  on  the  right  is  a  worn  sample  as  seen 
under  a  scanning  electron  microscope.  The  worn  surface  looks  fairly  smooth  except 
where  there  is  teardrop-shaped  plowing  damage. 

Wear  at  high  speeds  is  a  compilation  of  many  different  wear  mechanisms.  The 
sliding  contact  raises  the  temperature  of  the  metal  to  the  point  where  it  will  literally  melt 
away.  [5]  Lift  effects  on  the  sled  due  to  high  speeds  will  cause  the  slipper  to  “bounce” 
along  the  rail,  leading  to  inconsistent  surface  contact,  further  complicating  the  prediction 
of  thermal  and  frictional  properties  of  the  interacting  metal  surfaces.  Contact  interactions 
such  as  adhesion  and  abrasion  are  also  contributors  to  the  overall  wear.  An  example  of 
Adhesion  and  abrasion  can  be  seen  in  Figure  5. 
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Adhesion 


Abrosion 


Figure  5:  Adhesion  and  Abrasion  Wear  [6] 

Adhesion  occurs  when  two  surfaces  come  into  contact  and  stick  together.  If  the 
adhesion  is  strong  enough,  the  some  of  the  softer  material  is  removed  and  stays  gripped 
to  the  harder  material.  Eventually  the  extra  material  falls  off,  causing  additional  wear 
damage  in  the  process.  Abrasion  occurs  when  asperities  collide  with  enough  force  to 
fracture  the  asperity  and  cause  material  to  sheer  off.  Melt  wear  occurs  when  the  sliding 
materials  generate  so  much  heat  that  the  melting  temperature  of  one  or  both  of  the  metals 
is  reached.  Material  drips  off,  alleviating  some  of  the  friction  between  the  two  metals,  but 
as  mentioned  earlier,  not  decreasing  the  wear  overall. 

All  of  these  factors  and  other  tribological  conditions  make  predicting  the  slipper 
wear  difficult.  For  practical  and  financial  reasons,  it  is  desirous  to  be  able  to  model  the 
phenomenon.  Major  wear  incidents  can  lead  to  heavy  damage  to  both  the  slipper  and  the 
rail,  potentially  disrupting  expensive  experiments  and  requiring  additional  track  repairs. 
Not  being  able  to  predict  wear  also  makes  design  of  similar  high  speed  sliding 
components  difficult,  thus  a  way  to  determine  the  material  interaction  is  valuable  not  only 
to  those  at  HHSTT,  but  also  anyone  wanting  to  try  similar  high  speed  experiments. 
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1.2  Objective 

The  wear  problem  in  question  is  a  complicated  one  with  many  factors 
contributing  to  the  overall  loss  of  material.  It  is  of  interest  to  be  able  to  simplify  this 
problem  and  create  a  basic  model  that  will  predict  wear  accurately.  This  research  focuses 
on  developing  a  simple  plane  strain  model  using  a  hydrocode  program  called  CTH  to 
represent  a  3D  wear  phenomenon.  The  results  are  carried  out  without  the  initial  input  of  a 
temperature  history,  and  will  not  be  exact.  The  aim  is  see  if  a  simple  plane  strain  analysis 
has  some  merit  for  evaluating  the  wear  effects  at  high  velocity. 

1.3  Prior  Work 

This  is  definitely  not  the  first  investigation  into  the  wear  of  sliding  materials. 
There  have  been  many  attempts  at  wear  characterization  for  low  speeds  and  high  speed 
starting  as  early  as  the  1950s.  [7]  Research  into  the  wear  problem  at  HHSTT  began 
around  1960,  and  AFIT  has  picked  up  the  general  problem  around  1998.  [8,  9]  Much 
progress  has  been  made  in  this  area,  and  has  been  built  into  this  project.  The  following  is 
a  brief  review  of  the  research  and  conclusions  that  have  been  made  which  provide  the 
foundation  for  this  work. 

1.3.1  Wear  Research 

There  are  many  different  aspects  of  wear  to  investigate.  Concern  about  wear 
problems  first  arose  in  the  US  Army  when  wear  from  projectiles  in  cannons  was  taking  a 
toll  on  weapons  performance.  Much  progress  was  made  when  investigating  the  problem, 
though  research  on  this  subject  was  not  released  until  1975.  [10]  Early  work  proposed 
that  wear  was  primarily  a  function  of  load  and  mostly  independent  of  contact  area.  [7] 
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Investigation  began  in  the  1950s  into  the  effect  of  thermal  environments  in  wearing 
materials  and  found  that  if  the  temperature  of  the  material  rose  to  the  melting 
temperature,  an  interfluous  layer  forms  between  the  sliding  materials  and  is  removed  as 
the  sliding  progresses.  [5,  11,  12,  13]  This  causes  the  coefficient  of  friction  to  lesson 
between  the  materials,  but  does  not  necessarily  result  in  the  reduction  of  the  overall  wear 
rate. 

In  1979,  Challen  and  Oxley  introduced  the  idea  that  asperities  (microscopic 
surface  variations  on  material  faces)  can  be  used  to  create  models  for  wear,  noting  that 
the  contact  area  between  surfaces  is  really  just  the  area  of  asperities  interfacing  between 
the  materials.  [14]This  is  particularly  of  interest  to  this  research  since  asperities  are  used 
in  the  developed  model. 

It  was  becoming  apparent  that  there  were  different  factors  that  caused  different 
types  of  wear  to  dominate  wear  scenarios.  In  1987,  Lim  and  Ashby  began  developing 
“wear  maps”  that  showed  what  kind  of  wear  was  dominating  based  on  velocity  and  load. 
The  maps  predicted  when  particular  types  were  occurring,  such  as  delaminated  wear 
(mild  and  severe),  mild  oxidation  wear,  severe  oxidation  wear,  and  melt  wear.  [15] 
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Figure  6:  Lim  and  Ashby's  Wear  Map  for  Steel  [15] 

Further  progress  was  made  with  asperities  in  1982  with  research  done  by  Barber 

and  Bauer,  which  showed  that  at  certain  velocities,  asperity  behavior  changed  from  a 

quasi-steady  mechanical  behavior  at  low  speeds,  to  higher  speeds  when  inertia  effects 

began  to  dominate.  [16] 

For  as  long  as  wear  events  had  been  studied,  many  analytical  models  for 
predicting  wear  have  been  produced.  But  as  Meng  and  Ludema  found  in  their  study  in 
1995,  most  of  those  models  were  not  practical  for  widespread  use.  [17]  However,  with 
the  rise  of  numerical  methods  and  computational  modeling,  the  problem  becomes  much 
more  approachable.  Constitutive  models  that  are  able  to  track  strain  rate  were  introduced 
and  Finite  Element  Analysis  (FEA)  models  began  to  come  onto  the  scene  in  the  1990s. 
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[18,  19]  Additionally,  new  experimental  methods  have  allowed  for  better  characterization 
of  materials  at  high  speed  impacts.  [20] 

1.3.2  Holloman  Gouging  Research 

The  wear  problem  at  Holloman  AFB  has  its  own  history  as  well.  In  1960, 

Wolfson  began  research  into  high  speed  wear  utilizing  a  set  up  where  a  sled  was 
accelerated  down  a  track  to  a  specific  velocity.  [21]  Once  the  target  velocity  was  reached, 
a  pin  was  placed  in  contact  with  the  rail  for  a  specific  distance.  Wolfson  was  able  to 
calculate  wear  rates  by  measuring  and  weighing  the  worn  slider  and  comparing  it  to  its 
initial  state  prior  to  the  run.  Some  of  his  results  are  applicable  to  our  studies,  and  are  used 
as  comparative  data  in  this  research. 

The  gouging  phenomenon  was  first  observed  and  recorded  by  Graff  and  Dettloff 
in  1969,  when  experiments  progressed  to  speeds  around  1800m/s.  Teardrop-shaped 
gouges  occurred  randomly  along  the  rail  and  on  the  rocket  shoe.  [8]  The  gouges  were  2- 
4in.  long,  lin.  wide,  and  l/16in.  deep.  Graff  and  Dettloff  were  able  to  successfully 
ballistically  reproduce  the  phenomenon  in  the  lab.  Similar  gouging  wear  occurred  at 
Sandia  National  Labs  in  1973  as  observed  by  Gerstle,  who  was  able  to  perform  a 
metallurgical  examination  of  the  damaged  metal  and  identified  that  there  was  melt  wear, 
catastrophic  thermoplastic  shearing,  and  the  damage  affects  which  went  deeper  into  the 
metal  than  just  at  the  surface.  [22] 

In  2000,  HHSTT  began  the  development  of  a  new  track  designed  to  handle 
experiments  that  would  break  land  speed  records.  [23]  With  this  goal  in  mind, 
minimizing  wear  and  fully  understanding  the  dynamic  environment  of  the  track  became  a 
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high  priority.  Software  was  developed  called  the  Dynamic  Analysis  and  Design  System 
(DADS)  which  allowed  dynamic  time  dependent  characterization  of  the  experimental 
results  to  be  recorded.  This  data  has  been  instrumental  for  the  development  and  validation 
of  models  at  the  Air  Force  Institute  of  Technology  (AFIT). 

1.3.3  AFIT  Research 

There  has  been  an  ongoing  effort  at  AFIT  to  characterize  high  speed  wear  events, 
particularly  focusing  on  those  at  HHSTT.  Work  began  as  early  at  1998  with  Schmitz, 
where  the  gouging  problem  was  first  investigated.  [9]  Since  the  phenomenon  was  a  high 
energy  impact  problem,  whose  results  had  been  reproduced  with  ballistics,  it  was 
determined  that  the  hydrocode,  CTH,  developed  by  Sandia  specifically  for  ballistic 
events,  would  be  an  appropriate  program  to  model  the  high  energy  scenario.  Laird  in 
2002  continued  this  research  and  developed  a  model  in  CTH  that  had  a  flyer  (slipper) 
colliding  with  a  target  (rail).  [24]  He  experimented  with  this  model  using  various 
temperature  profiles,  horizontal  and  vertical  velocities,  collisions  with  asperities, 
collisions  without  asperities,  and  considered  a  gap  between  the  rail  and  the  slipper. 
Szmerekovsky  in  2006  continued  the  investigation  while  adding  coatings  to  the  rail  to  see 
how  they  would  affect  the  material  surface  temperatures.  [25]  Also  in  2006,  Cinnamon 
had  a  chance  to  do  a  metallographic  analysis  on  materials  from  previous  tests  at  the 
HHSTT,  and  developed  improved  material  data  for  the  model.  [4]  In  2007,  Cameron 
attempted  to  separate  the  overall  wear  into  two  different  categories-  melt  and  mechanical. 
He  used  the  CTH  code  to  calculate  each  event  separately  and  was  able  to  compare  his 
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results  to  a  2003  experiment  data.  His  results  were  found  to  be  reasonably  equivalent. 

[26] 

In  2008,  Chmiel  moved  away  from  the  CTH  code  and  embarked  on  developing  a 
finite  element  code  in  ABAQUS  to  model  wear.  [27]  He  attempted  two  different  models: 
one  macroscopic  and  one  microscopic.  The  macroscopic  code  experienced  technical 
difficulties,  but  the  microscopic  code,  employing  specific  material  characteristics  and 
failure  criteria,  produced  favorable  results.  This  work  was  expanded  on  by  Hale  in  2009, 
which  took  the  plane  strain  ABAQUS  model  and  was  able  to  use  it  to  compare  to  a  2008 
HHSTT  experiment.  [3] 

The  CTH  model  was  not  left  behind  with  the  rise  of  the  FEA  model  as  it  was  still 
computationally  superior  to  FEA  codes  at  higher  speeds.  Using  Lim  and  Ashby’s 
research,  Meador  continued  using  the  CTH  model  and  attempted  to  determine  if  the  heat 
distribution  coefficient  between  the  slipper  and  the  rail  could  be  approximated  at  50% 
(0.5).  He  found  that  since  the  slipper  contact  with  the  rail  was  not  consistent,  0.5  was  a 
poor  assumption  and  more  reasonable  to  assume  a  coefficient  of  0.12-0.14.  [28]  In  2011, 
Huber  continued  development  of  the  CTH  code  by  experimenting  with  different  failure 
criteria  to  see  how  they  would  affect  the  predicted  wear  results.  [29]  Work  is  currently 
ongoing  with  Buentello  continuing  the  work  on  the  ABAQUS  FEA  model. 

1.4  Prior  Work  Summary 

Research  into  wear  has  revealed  that  it  is  a  complicated  mechanism  with  several 
different  types  occurring  simultaneously,  and  the  wear  behavior  varies  based  on  speed.  At 
speeds  over  1828.8  m/s,  it  has  been  found  that  teardrop  gouging  develops  during  testing 
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at  HHSTT,  and  these  same  results  have  been  reproduced  in  the  lab  using  ballistic 
impacts.  At  AFIT,  finite  element  models  (both  2D  and  3D)  have  been  developed  in 
ABAQUS.  Additionally,  a  2D  CTH  code  using  various  failure  criteria  has  been 
developed.  The  models  are  in  the  process  of  being  refined,  and  are  validated  using 
experimental  data  from  HHSTT  experiments. 


II.  Theory 


2.1  Introduction 

On  a  microscopic  scale,  surfaces  are  not  smooth.  They  are  actually  composed  of 
uneven  jagged  profiles,  contact  between  the  interfaces  of  the  two  surfaces  only  occurs  at 
the  tips  of  which  touch  the  other  surface.  Figure  7  shows  a  representation  of  the  uneven 
surface  profile. 


20 


20 
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Figure  7:  Surface  Profile  of  a  Metal  [30] 

These  contact  points  can  be  modeled  as  small  bumps,  or  “asperities.”  At  high 
speeds,  the  collisions  between  asperities  are  actually  impact  phenomena,  meaning  that 
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there  is  the  possibility  of  internal  shock  wave  propagation.  Inertia  effects  and  stress  wave 
propagation  become  significant  factors  at  play.  [31,  32],  thus  the  employment  of  classic 
deformable  body  mechanics  which  use  Hooke’s  Law  to  predict  stress,  strain,  and 
deformation  are  not  enough  to  explain  what  is  happening  with  stress  propagation  inside 
the  materials,  since  materials  will  behave  differently  at  high  speeds  than  at  low  speeds, 
and  a  different  approach  to  modeling  the  event  must  be  used. 

2.2  Theory  Overview 

This  problem’s  solution  revolves  around  understanding  how  stress  waves  will 
travel  through  materials,  along  with  trying  to  satisfy  the  three  conservation  equations:  the 
conservation  of  mass,  the  conservation  of  momentum,  and  the  conservation  of  energy. 
High  speed  impact  complicated  things,  since  materials  will  become  more  compressible 
and  due  to  discontinuities  in  pressure  in  the  wave  profile,  the  thermodynamic  behavior  of 
the  material  is  not  standard.  [33]  Thus,  more  relations  are  needed  to  complete  the  solution 
set.  An  equation  of  state  is  introduced  to  take  the  unique  density  and  thermodynamic 
behavior  of  the  materials  into  account.  A  constitutive  equation  rounds  out  the  solution 
and  defines  how  the  stress  developed  over  time  taking  into  account  strain,  energy,  and 
other  factors  depending  on  the  model  that  is  used.  Additionally,  predicting  dynamic 
failure  is  a  challenge,  so  failure  criteria  are  developed  in  an  attempt  to  characterize  it. 

While  the  conservation  equations  are  standard  for  all  solutions,  models  can  be 
varied  based  on  the  selection  of  the  equations  of  state  that  represent  the  materials,  the 
constitutive  models  that  develop  stress,  and  for  numerical  modeling,  on  how  the  program 
determines  failure  of  a  material  during  simulations.  Thus,  there  numerous  possible 
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models  using  various  combinations  of  these  components.  These  models  can  be  validated 
by  comparing  them  to  experimental  data,  and  it  is  attempted  to  do  so  in  this  research. 

The  conservation  equations,  an  equation  of  state,  and  a  constitutive  equation,  and 
a  specific  failure  criterion  are  the  basic  tools  that  are  used  to  solve  this  problem.  First,  the 
basics  of  wave  propagation  are  explored  to  better  understand  why  the  thermal  properties 
are  discontinuous.  Then  the  specifics  of  each  of  the  relationships  are  further  detailed  in 
the  following  sections. 

2.3  Wave  Propagation 

Meyer’s  Dynamic  Behavior  of  Materials  gives  an  excellent  overview  of  wave 
propagation.  Much  of  Section  2.3  is  a  summary  of  chapters  1-4  of  his  book.  [34] 

Technically  speaking,  everything  is  moving  in  waves.  Atoms  are  constantly 
vibrating,  but  for  most  situations,  the  vibrations  are  too  small  to  notice.  They  exist  in  a 
static  or  quasistatic  equilibrium.  Application  of  forces  can  be  described  as  vibrations  as 
well.  Static  deformation  revolves  around  balancing  forces  in  equilibrium  without  much 
thought  to  time  and  how  stress  moves  through  the  material,  because  the  application  of 
force  can  be  considered  instantaneous.  This  is  not  true  however  in  dynamic  deformation 
where  bodies  are  subjected  to  rapidly  changing  loads.  In  these  cases,  one  must  be  able  to 
track  the  velocity  of  a  stress  wave  moving  through  a  body. 

To  begin  to  understand  how  stress  waves  travels,  it  is  helpful  to  understand  how  a 
simple  uniaxial  sinusoidal  wave  operates.  A  low,  constant  amplitude  wave  advances 
according  to  equation  2. 1 : 
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v  —  Af 


2.1 


Where  v  is  the  velocity,  A  is  the  wave  length,  and  f  is  the  frequency. 

The  vertical  position  of  the  wave  can  be  defined  by  equation: 

y  —  To  sin  [27 rf(t  -  £)]  2.2 


Where  y  is  the  displacement,  yo  is  the  wave  amplitude,  t  is  the  time,  and  x  is  the 
position. 

Now,  assume  that  the  wave  can  be  represented  by  a  string  with  a  mass  m  and  a 
line  tension  T.  To  get  a  wave  equation  in  the  differential  form,  one  carries  out  an 
equilibrium  analysis  of  the  ends  of  a  wave  and  applies  Newton’s  Law  (assuming  all 
motion  is  in  the  Y  direction): 


2.3 


p  ds 


d2y 
dt 2 


2.4 


d2y 

Where  (p  ds)  is  the  mass  per  unit  length,  -^-7  is  the  acceleration. 

Summing  the  forces  at  the  ends,  accounting  for  the  tension  inside  of  the  string, 
and  assuming  small  oscillations  (cosine  theta  =1),  brings  us  to  the  differential  wave 
equation. 


d2y  _  2  d2y 

P  dt2  dx2 


2.5 


Where  T  is  the  tension  in  the  string. 

Differentiate  position  equation  twice  and  substitute  it  into  the  differential  wave 
equation  gives  the  traveling  wave  velocity  relation: 
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2.6 


i 


In  the  next  section,  it  will  be  shown  that  the  tension  is  replaced  by  the  modulus  of 
elasticity  for  elastic  waves.  Substituting  equation  2.6  into  2.5,  the  differential  wave 
equation  for  a  horizontal  traveling  wave  with  a  horizontal  displacement  becomes: 

2.7 

dt2  dx2 

This  equation  applies  when  any  shape  of  disturbance  propagating  at  velocity  v  in 
a  medium  and  remaining  unchanged.  It  is  an  analogous  representation  to  the  development 
of  traveling  elastic  waves.  Each  traveling  wave  in  a  body  is  going  to  do  so  at  its  own 
characteristic  velocity.  Velocities  in  traveling  waves  are  a  function  of  the  stress  response 
of  the  material  to  strain  and  density,  and  so  will  be  different  based  on  whether  the  wave  is 
elastic,  plastic,  or  includes  shocks. 

2.3.1  Elastic  Waves 

For  elastic  waves,  the  stress-strain  relations  are  based  on  elastic  constants.  When 
the  dynamic  process  (application  of  force)  is  slow,  a  body  can  be  considered  in  static 
equilibrium,  however  this  is  not  technically  the  case.  Internal  stresses  are  not  immediately 
transferred  though  materials,  they  travel  atom-to-atom  via  waves.  Essentially  all  atomic 
movement  can  be  traced  back  to  the  conservation  of  momentum  equation.  The  next 
section  explains  how  elastic  stress  travels  in  a  material. 

For  elastic  waves  in  a  bar,  deformation  is  elastic  and  follows  Hooke’s  Law: 

a  —  Ee  2.8 

Where  a  is  the  stress,  E  is  the  elastic  modulus,  and  e  is  strain. 
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Horizontal  strain  can  be  defined  in  terms  of  displacement: 


£ 


du 

dx 


2.9 


Where  u  is  displacement. 
An  elastic  wave  travels  as 


d2u  E d2u 
dt2  p  dx2 
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Comparing  this  to  equation  2.1 1,  the  velocity  of  wave  is  defined  as: 


C0  = 


2.11 


There  are  several  types  of  elastic  waves  that  may  occur.  Longitudinal  (or 
irrotational  or  dilatational)  waves  occur  in  infinite  and  semi-infinite  media,  and  travel  in 
the  direction  of  the  initial  contact.  Distortional  (or  shear  or  transverse  or  equivolumal) 
waves  occur  in  solids  perpendicular  to  the  initial  contact  direction,  and  there  is  no  change 
in  density  or  longitudinal  strains  Sn,  822,  and  S33.  Surface  (or  Rayleigh)  waves  move  up 
down  back  and  forth  along  the  surface.  Interfacial  (or  Stanley)  waves  occur  when  two 
semi-infinite  materials  with  different  properties  are  in  contact  with  each  other.  Finally 
bending  waves,  which  occur  in  bars  and  plates,  in  2D  scenarios  are  the  waves  caused  by 
bending  stresses.  Each  kind  of  wave  travels  at  a  different  velocity:  surface  waves  are 
slow,  while  longitudinal  waves  are  fast. 

At  boundaries,  elastic  waves  divide  into  reflection  and  refraction  when  they  come 
into  contact  with  mediums  with  different  sonic/sound  impedances.  The  specific  angles  of 
the  reflected  and  refracted  waves  can  be  found  using  Lenz’s  law.  These  waves  are  easiest 
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to  deal  with  when  the  impact  is  normal  to  the  surface  in  question.  The  total  initial 
stress/particle  velocity  is  equal  to  the  sum  of  the  transmitted  and  reflected  velocities.  The 
impedance  of  materials  determines  the  amplitudes  of  the  resulting  transmitted  and 
reflected  waves.  The  impedance  is  based  on  the  material  properties  of  the  material,  and 
the  relationship  is  seen  in  equation  2.12. 

Impedance  —  pC0  2.12 

Where  p  is  the  density  and  C0  is  the  elastic  wave  velocity. 

Incident  wave 


Figure  8:  Materials  with  Differing  Impedances 

For  example,  consider  two  materials  are  subjected  to  an  incident  wave  as  shown 
in  Figure  8.  When  impedance  of  blue  material  is  greater  than  the  red  material,  a  pulse  of 
the  same  sign  is  reflected.  When  impedance  of  the  red  material  is  greater  than  blue 
material ,  a  pulse  of  the  opposite  sign  is  reflected. 

On  a  free  surface,  the  elastic  modulus  (E)  is  equal  to  0.  On  a  rigid  body,  E  is  equal 
to  oo.When  a  compressive  wave  encounters  a  free  surface,  it  reflects  back  as  a  tensile 
wave.  The  stress  sign  changes  but  the  velocity  stays  the  same.  When  a  compressive  wave 
encounters  a  rigid  body,  a  compressive  wave  is  reflected  back  (the  stress  sign  stays  the 
same  but  the  velocity  changes).  If  the  particle  velocity  and  the  wave  direction  have  the 
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same  direction,  the  stress  is  compressive.  If  they  have  opposite  senses,  the  stress  is 
tensile. 

2.3.2  Plastic  Waves 

For  plastic  waves,  the  stress-strain  relations  are  based  on  work  hardening  and 
strain  rate.  Plastic  waves  occur  when  a  stress  in  a  material  exceeds  its  elastic  limit.  This 
happens  in  both  static  and  dynamic  deformations.  For  elastic  waves,  the  stress  was 
entirely  elastic  and  could  be  based  on  a  linear  equation.  This  is  not  necessarily  the  same 
for  plastic  stress.  Figure  9  below  shows  4  different  stress-strain  curves  and  how  they 
change  depending  on  how  plastic  flow  stress  develops. 
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Figure  9:  Stress  Strain  curves  for  ductile  materials  with  (a)  bilinear  elastoplastic;  (b)  power  law  work 
hardening;  (c)  strain  rate  dependent  flow  stress;  (d)  strain  rate  history  dependent  flow  stress 

In  Figure  9  (a),  one  assumes  that  both  the  plastic  and  the  elastic  portions  of  the 
stress-strain  curve  are  assumed  to  be  linear  relations.  In  Figure  9  (b),  the  plastic  portion  is 
assumed  to  develop  based  on  a  power  law,  where  the  power  exponent  is  the  work 
hardening.  In  Figure  9  (c),  the  stress  is  increasing  as  the  applied  strain  rate  increases. 
Finally,  in  Figure  9  (d),  the  stress  is  increasing  with  the  increasing  strain  rate,  but  is  also 
dependent  on  when  the  strain  rate  is  applied.  One  general  stress-strain  relation  for  a 
changing  strain  rate  is  represented  as: 

a  =  a0  +  k£nem  2.13 
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Where  er0  is  the  material  yield  stress,  n  is  work  hardening,  and  m  is  the  strain  rate 
sensitivity. 

Temperature,  strain,  and  strain  rate  all  affect  mechanical  material  behavior  and 
how  the  flow  stress  develops. 

To  get  an  idea  of  how  the  velocity  of  a  plastic  wave  differs  from  that  of  an  elastic 
wave,  let’s  look  at  the  development  of  stress  in  a  plastic  wave.  For  time  independent 
waves,  the  speed  that  they  travel  at  can  be  found  to  be: 


2.14 


20 


There  are  three  regions  -  where  the  plastic  wave  is  dominating,  where  the  elastic 
wave  is  dominating,  and  then  when  the  wave  comes  to  a  complete  stop.  Note  how  CO  is 
faster  than  Cl,  but  Cl  incurs  more  strain. 

To  find  the  velocity  at  x  =  0  (right  when  the  plastic  wave  initiates),  one  arrives  at 
this  integral  in  equation  2.15: 

der0 

3Lde  2A5 

Po 

Where  Vi  is  the  velocity  at  x=0  and  Si  is  the  maximum  plastic  strain. 

This  equation  must  be  solved  numerically,  and  it  is  impossible  to  know  the 
velocity  at  this  point  without  knowing  the  relationship  between  stress  and  strain.  This 
relationship  is  known  as  a  constitutive  model.  There  have  been  many  models  developed 
and  they  have  been  divided  into  several  different  categories:  linear  viscoelastic, 
infinitesimal  elastic-perfectly  plastic,  nonlinear  viscoelastic,  and  thermoelastic. 

2.3.3  Shock  Waves 

Shock  waves  that  occur  in  a  material  are  largely  due  to  adiabatic  compressibility. 
Shear  stresses  are  neglected  since  the  event  is  at  a  high  pressure.  When  the  amplitude  of 
the  stress  wave  is  greatly  exceeding  the  flow  strength  of  a  material,  a  disturbance  front 
will  “steepen  up”  as  it  travels  through  the  material.  This  is  because  as  the  strain  increases, 
higher  amplitude  regions  of  the  front  travel  faster  than  the  lower  regions.  This  leads  to  a 
shock  wave  -  “a  discontinuity  in  pressure,  temperature  (or  internal  energy),  and  density.” 
[34]  This  is  illustrated  in  Figure  1 1 . 
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Distance  (or  time)  — ► 


Figure  11:  Increasing  Pressure  with  Time  [32] 

The  wave  is  shown  traveling  as  one  would  expect,  with  the  bottom  part  of  the 

wave  with  lower  pressure,  representing  where  elastic  stress  would  be  occurring,  is 
moving  faster  than  the  upper  part  of  the  wave  where  presumably  plastic  deformation 
would  occur.  As  the  wave  progresses,  the  upper  part  of  the  wave  speeds  up  and  catches 
up  with  the  lower  part  of  the  wave.  Now  the  pressure  wave  no  longer  looks  like  a  wave, 
but  rather  a  straight  vertical  line  causing  a  discontinuity  in  the  pressure  data,  causing  it  to 
no  longer  be  a  continuous  function.  Because  of  this,  additional  shock  relations  must  be 
employed  to  get  us  across  the  jump.  This  is  possible  by  re-working  the  conservation 
equations  to  track  the  different  particle  velocities  and  employing  an  equation  of  state  to 
create  thermodynamic  relationships  of  the  materials  at  high  speeds. 

2.4  Conservation  Equations 

To  solve  this  problem,  one  approaches  wielding  the  basic  conservation  equations: 
conservation  of  mass,  conservation  of  momentum,  and  conservation  of  energy,  as  shown 
below: 
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Conservation  of  mass  shows  that  material  is  neither  created  nor  destroyed. 


I  pdV  —  const  2.16 

'v 

Where  p  is  the  material  density  and  V  is  the  material  volume. 

The  conservation  of  momentum  equation  shows  that  the  net  force  in  a  system  is 
equal  to  the  rate  of  change  of  momentum. 


F  = 


2.17 


Where  F  is  force,  m  is  mass,  v  is  velocity,  and  t  is  time 
Conservation  of  Energy  states  that  the  total  overall  energy  in  a  system  is 
conserved. 


'YdEi='Yd\pvf ='YdEf+'Yd\pvf+w  2as 

Where  E  is  the  internal  energy,  -pv  is  the  kinetic  energy,  W  is  the  work  done  on 
the  system,  and  i  and  f  are  indices. 

These  equations  track  energy,  pressure,  density,  and  velocity  inside  of  the 
materials.  Due  to  the  compressibility  caused  at  high  speeds,  the  density  may  not  remain 
constant  throughout  the  event.  Additionally,  there  are  the  discontinuities  caused  by 
possible  shocks  to  take  into  account.  Consider  the  one-dimensional  case  illustrated  below 
in  figure  12: 
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Figure  12:  Compression  of  a  Medium  with  a  Generated  Shock  Wave  [34] 

Figure  12  (a)  shows  the  initial  state  of  shock  impact  at  time  (t)  =  0,  where  the 

material  in  question  is  at  rest.  In  figure  12  (b),  the  plunger  moves  forward,  compresses 

the  material,  and  sends  a  shockwave  forward  ahead  of  the  initial  impact  site. 

Rankine  and  Hugoniot  went  about  rearranging  the  conservation  equations  for  a 

ID  case,  and  are  as  displayed  below: 


Pq^S  pijj  S  Up') 

2.19 

P-P0=  Po(UsUp) 

2.20 

,  ,  1  / 1  1\ 

2.21 

Where:  pQ  is  the  initial  density,  Us  is  the  shock  velocity,  p  is  the  density,  Up  is  the 
particle  (impact)  velocity,  P  is  the  pressure,  P0  is  the  initial  pressure,  E  is  the  internal 
energy,  and  E0  is  the  initial  internal  energy. 
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The  three-dimensional  equations  are  more  complicated,  but  based  on  the  same 
principles. 

Note  from  equations  2.19-2.21  that  there  are  five  variables  of  interest:  density, 
shock  velocity,  particle  velocity,  pressure,  and  energy.  Since  there  are  more  unknowns 
than  equations  to  solve  them,  additional  relationships  are  needed.  The  first  three 
relationships  are  provided  by  the  conservation  equations,  and  the  fourth  relationship  is 
provided  by  an  equation  of  state.  With  the  addition  of  an  equation  of  state,  the  ability  is 
gained  through  experimentation  to  relate  four  of  the  unknowns  in  terms  of  one  unknown. 

2.5  Equation  of  State 

The  equation  of  state  helps  relate  the  thermodynamic  properties  of  a  material  to 
its  density  and  pressure.  [31]  They  help  deal  with  pressures  that  are  higher  than  material 
strength,  compressibility  effects,  and  shock  heating.  Equations  of  state  are  dealing  with 
the  hydrodynamic  (volumetric)  portion  of  stress,  while  constitutive  equations  deal  with 
the  deviatoric  portion  of  stress.  They  are  of  the  general  form: 

P  =  P(p,E )  2.22 

Where  P  is  the  pressure,  p  is  the  density  and  E  is  the  internal  energy. 

Equation  of  states  can  be  either  tabular  or  analytic. 

2.5.1  SESAME  tables 

The  equation  of  state  used  in  this  study  is  a  set  of  tables  called  the  SESAME 
tables.  It  consists  of  thermodynamic  properties  such  as  pressure,  density,  and  energy 
experimentally  determined  by  Los  Alamos  National  Laboratory.  [35]  Additionally, 
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vaporization,  melt,  and  shear  data  is  available,  as  well  as  shock  data.  This  is  notable 
because  it  is  inclusive  of  any  phase  changes  that  are  occurring. 

2.5.2  Mie-Griineisen 

The  Mie-Griineisen  Equation  of  state  is  an  example  of  an  analytic  EOS.  It  relates 
the  materials  pressure  and  internal  energy  to  its  pressure  and  temperature  with  the 
following  equation: 


Where: 


V  = 


PqCqV  ( 
(1  -  srj)2  V 


2.23 


1  p0 
n  =  1 - 

P 

2.23 

Us  =  c0  +  sU0 

2.24 

T  is  the  Griineisen  parameter,  Co  is  the  reference  speed  of  sound,  E  is  the  internal  energy, 
s  is  the  parameter  that  relates  shock  and  particle  velocity,  p  is  the  density,  and  p  is  the 
pressure. 

It  allows  for  a  finer  control  of  temperature  inputs,  but  isn’t  set  up  to  handle  phase 
changes. 


2.6  Constitutive  Equation 

Stress  can  be  divided  into  two  components  -  a  volumetric  portion  generating 
stresses  from  volumetric  changes  and  a  deviatoric  portion  generating  stresses  from 
changes  in  shape  of  the  body.  [36] 

Oij  -  -PSij  +  stJ  2.25 
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Where  cri;  is  he  total  stress,  —  is  the  volumetric  stress,  and  sLj  is  the 


deviatoric  stress  component. 

The  volumetric  stresses  can  be  found  using  data  from  the  equation  of  state  and 
hydrostatic  stress  (the  sum  of  the  principal  stresses),  but  an  additional  relation  is  needed 
to  find  the  deviatoric  stress.  This  relation  comes  from  a  constitutive  equation,  which 
calculates  the  dynamic  yield  stress  (or  flow  stress)  in  terms  of  strain  rates  and 
temperature.  [37].  Note  that  the  dynamic  yield  strength  is  not  the  same  as  the  static  yield 
strength.  Additionally,  the  dynamic  flow  stress  is  not  representative  of  the  entire  stress- 
strain  curve,  but  only  the  portion  where  plastic  deformation  has  been  initiated. 
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This  relationship  is  at  a 
particular  value  of  strain. 
Increasing  strain  would  shift  the 
curves  higher,  as  increasing 
strain  yields  higher  material 
strength  (i.e.  strain-hardening). 
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Figure  13:  General  Constitutive  Model  Behavior  [37] 

Figure  13  shows  the  general  behavior  of  constitutive  models  and  how  they  change 

based  on  increasing  strain  rate  and  temperature.  Only  strain  rate  and  temperature  changes 
are  accounted  for  in  this  example,  but  models  can  take  into  account  other  factors,  such  as 


27 


material  properties  and  work  hardening.  In  addition  to  modeling  the  deviatoric  stress, 
constitutive  models  help  determine  when  a  material  is  behaving  elastically  or  plastically. 
An  equivalent  stress  is  calculated  using  the  second  invariant  of  the  deviatoric  stress 
tensor: 


°EQ  ~ 


J 


3  (Si  jsi  j ) 


2.26 


When  the  equivalent  stress  is  less  than  the  flow  stress,  the  material  is  assumed  to 
be  in  its  elastic  deformation  stage.  When  the  equivalent  stress  is  greater  than  the  flow 
stress,  the  material  is  assumed  to  be  in  a  state  of  plastic  deformation.  [36] 

The  constitutive  model  used  in  this  research  is  the  Johnson  and  Cook  Viscoplastic 
Flow  equation.  “Viscoplastic”  indicates  that  the  material  behavior  is  rate  dependent  [38]. 
The  equation  is  shown  below  as:  [39] 

<7  =  (A  + Ben)(l  + C\ne*)(l-T*m)  2.27 

Where: 


A  =  initial  yield  strength 
B  =  hardening  modulus 
C  =  strain  rate  dependent  coefficient 
m  =  thermal  softening  coefficient 
n  =  strain  hardening  coefficient 
£  =  equivalent  plastic  strain 

s*  =  the  equivalent  plastic  strain  rate/effective  plastic  strain  rate 
T*  =  homologous  temperature  =  ^  Troom) 

T  melt~Troom 

T  =  Absolute  temperature 

The  equation  is  divided  into  three  bracketed  components  -  the  first  set  of  brackets 
takes  into  account  the  strain  when  T*  =  0  and  s*  =  1,  the  second  set  of  brackets  takes  into 
account  the  strain  rate  effects,  and  the  third  set  of  brackets  takes  thermal  effects  into 
account. 
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Coefficients  for  this  model  come  from  experimental  testing  using  both  the  Split  - 
Hopkinson  bar  test  and  flyer  plate  experiments  which  allow  material  properties  to  be 
obtained  within  dynamic  conditions. 

2.6.1  Split  Hopkins  on  Bar  Test 

The  Split-Bar  Hopkinson  test  is  carried  out  by  placing  a  material  sample  in 
between  two  bars  -  the  incident  bar,  and  the  transmission  bar.  A  schematic  of  the 
experiment  is  shown  below  in  Figure  14. 


Figure  14:  Split-Hopkinson  Bar  Test  Apparatus  Schematic  [37] 

The  incident  bar  is  struck  by  the  striker  bar  causing  an  elastic  pulse  wave  to  travel 

through  the  specimen  and  the  transmission  bar.  Strain  gauges  set  up  on  the  incident  bar 
and  the  transmission  bars  are  able  to  measure  the  strain  waves  that  are  transmitted  and 
reflected  [40].  Cinnamon  in  his  dissertation  was  able  to  perform  this  test,  and  his 
experimental  set  up  is  shown  in  Figure  15. 
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Figure  15:  Split  Hopkinson  Bar  Apparatus  Test  Section  at  the  University  of  Dayton  [37] 


Using  his  results,  he  was  able  to  determine  the  physical  properties  and  material 
constants  used  in  the  Johnson  and  Cook  Constitutive  Model.  His  results  are  shown  in 
Table  1. 

Table  1:  Cinnamon's  Experimental  Results  in  the  Physical  Properties  and  Model  Constants  from  the 

Split-Hopkinson  Bar  Test  [37] 


Property /Const  ant 

1080  Steel 

VascoMax  300 

E  (GPa) 

202.8 

180.7 

V 

0.27 

0.283 

Tmclt  (K) 

1670 

1685 

p  ( kg/m 3) 

7800 

8000 

JC:  A  (GPa) 

0.525 

2.17 

JC:  D  (GPa) 

3.59 

0.124 

JC:  C 

0.029 

0.0046 

JC:  m 

0.7525 

0.95 

JC:  n 

0.6677 

0.3737 
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2.6.2  Flyer  Plate  Experiment 

The  Split-Hopkinson  Bar  tests  create  strain  rates  on  the  order  of  10  /second,  but 
higher  strain  rates  are  required  to  fully  characterize  the  constitutive  model  (on  the  order 
of  104/sec-105/sec).  Flyer  plate  impact  tests  reach  these  high  strain  rates,  and  allow  for  the 
stress  to  be  measured  with  respect  to  time.  In  flyer  plate  experiments,  it  is  the  goal  to 
strike  a  plate  uniaxially  with  a  flyer  plate.  A  schematic  of  the  experiment  is  shown  below 
in  Figure  16. 
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..  ■ 

Flyer 

Plate 

■ 

Gun  Barrel 
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Figure  16:  Flyer  Plate  Experiment  Schematic  [37] 

Stress  gauges  behind  the  plate  record  the  stress.  The  measurements  should  be  able 

to  pick  up  the  plastic  and  elastic  waves  occurring  and  determine  the  peak  stress  and 

Hugoniot  Elastic  Limit  as  seen  in  Figure  17.  [35] 
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Figure  17:  Idealized  Stress  vs.  Time  Plot  for  a  uniaxial  Planar  Impact  [37] 

Cinnamon  carried  out  these  experiments  in  2006.  Figure  18  shows  the  flyer  plate 

and  target  used  in  his  experiment. 


Figure  18:  Flyer  Plate  and  Target  used  in  University  of  Dayton  Experiments  [37] 

From  data  collected  from  these  tests,  all  the  constants  can  be  determined.  Results 

are  seen  in  Table  2. 
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Table  2:  Cinnamon's  Experimental  Results  in  the  Physical  Properties  and  Model  Constants  from  the 

Flyer  Plate  Tests  [37] 


Property  /  Constant 

1080  Steel 

VascoMax  300 

E  (GPa) 

202.8 

180.7 

u 

0.27 

0.283 

Tmcll  (K) 

1670 

1685 

p  (kg/m2) 

7800 

8000 

JC:  A  (GPa) 

0.7 

2.1 

JC:  D  (GPa) 

3.6 

0.124 

JC:  C 

0.017 

0.03 

JC:  m 

0.25 

0.8 

JC:  n 

0.6 

0.3737 

To  obtain  the  final  material  constants  to  be  used,  the  constants  were  adjusted  to 
better  fit  experimental  stress  data.  [37] 

2.6.3  Johnson  and  Cook  Constitutive  Model  Behavior 

A  plot  of  the  stress  developed  in  a  cell  using  ABAQUS  gives  some  insight  to  how 
the  stress  is  developing.  This  plot  is  seen  in  Figure  19. 
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Figure  19:  Johnson  and  Cook  Flow  Stress  Behavior  Generated  in  ABAQUS  [38] 

The  feature  of  interest  in  Figure  19  is  the  blue  line,  which  is  the  plot  of  the 

Johnson  and  Cook  flow  stress  curve.  The  beginning  of  the  curve  shows  the  elasticity 

portion  of  stress,  and  is  clearly  linear.  After  hitting  a  maximum  value,  the  plastic  flow 

stress  begins  to  develop.  As  the  temperature  and  equivalent  plastic  strain  increase,  the 

plastic  stress  decreases  until  it  hits  a  critical  failing  value.  The  stress  drops  off  from  there 

(In  CTH  the  stress  is  directly  set  to  0) 

The  fact  that  the  model  includes  the  increase  in  temperature  due  to  adiabatic 
heating  from  plastic  deformation  in  its  calculation  of  flow  stress,  takes  into  account  the 
thermal  softening  and  strain  rate  effects,  in  addition  to  having  the  material  constants  for 
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the  model  for  our  specific  materials,  makes  this  model  an  excellent  candidate  for  this 
research.  Additionally,  it  has  its  own  accompanying  failure  criterion  that  is  discussed  in 
the  next  section. 

2.7  Failure  Criteria 

The  purpose  of  failure  criteria  is  to  predict  what  conditions  will  cause  failure  of  a 
material.  There  are  multitudes  of  different  criteria  to  choose  from,  and  each  criterion  will 
give  different  results  for  when  failure  occurs  and  how  much  damage  accumulates.  Thus, 
the  selection  of  a  damage  criterion  is  a  very  influential  part  of  modeling  wear,  as  different 
criteria  will  predict  different  damage  levels,  and  thus  different  wear  levels.  The 
differences  in  wear  criteria  are  amplified  even  further  for  dynamic  fracture.  The  flow 
equations  are  more  complicated  than  in  quasistatic  fracture,  and  so  criteria  are  more 
difficult  to  reliably  develop  for  all  dynamic  impact  cases.  This  research  uses  two  failure 
criteria  developed  in  previous  research,  as  well  as  one  new  criterion  that  previously  could 
only  be  used  in  ABAQUS. 

2.7.1  Johnson  and  Cook  Failure  Criterion 

This  study  focuses  mainly  on  the  Johnson  and  Cook  Failure  Criterion,  which  is 
meant  to  be  used  in  conjunction  with  the  Johnson  and  Cook  constitutive  equation.  It 
dictates  when  a  material  has  failed  based  on  strain  rates.  [39]  It  can  be  written  as: 

f  de 

D  =  - v  2.28 

J  £Pf(T,p,£,£ ) 

Where: 

D=  Damage  coefficient 
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ds  =  changing  strain 

epf  =  equivalent  plastic  strain  at  fracture 

T  =  absolute  temperature 

p  =  pressure 

e  =  strain 

e  =  equivalent  deviatoric  strain  rate 


The  equivalent  plastic  strain  is  calculated  with  this  equation  wherein: 


£Pf  = 


D 1  +  D 2  exp  (-D3p)  [D4  +  ln(max(l,  e))][l  +  DST*] 


Where 

Y  =  yield  stress  (from  Johnson  and  Cook  model) 
p  =  pressure 

D1-D5  =  experimentally  determined  coefficients 


2.29 


The  equation  is  divided  into  three  components:  the  first  set  of  brackets  denotes 
pressure  dependence,  the  second  set  denoted  strain  rate  dependence,  and  the  third  denotes 
temperature  dependence.  [38]  The  coefficients  for  the  first  set  of  brackets  (D1-D3)  are 
found  experimentally  from  isothermal  torsion  tests,  while  D4  and  D5  for  the  strain  and 
temperature  effects  are  found  from  un-notched  tension  tests  and  Hopkinson  bar  tests  with 
varying  temperatures.  [39] 

The  damage  coefficient  (D)  ranges  in  value  from  0  to  1,  where  0  indicates  that  no 
failure  has  occurred  in  the  material,  and  1  indicates  total  failure  of  the  material.  Numbers 
in  between  indicate  partial  failure.  Once  the  material  has  failed,  its  tensile  strength  is  set 
to  0  and  it  loses  its  ability  to  support  shears.  [39]  This  is  a  unique  feature  of  the  Johnson 
and  Cook  failure  criterion  that  is  not  seen  in  other  criteria  used.  It  behaves  much  like  a 
fluid  and  is  essentially  removed  from  the  scenario. 


36 


2.7.2  Strain  at  Max  Stress  Failure  Criterion 


This  method  determines  failure  based  on  the  plastic  strain  at  maximum  stress  for  a 
given  strain  rate.  True  stress  verses  true  strain  was  plotted  by  Hale  for  several  strain  rates, 
as  shown  below  in  Figure  20. 


Figure  20:  True  Stress-Strain  Curves  for  Varying  Strain  Rates  using  the  Johnson  and  Cook 

Constitutive  Equation  [3] 

The  flow  stress  is  the  stress  calculated  using  the  Johnson  and  Cook  constitutive 
equation.  Each  curve  is  the  true  stress-  true  strain  plot  of  VascoMax300  at  various  strain 
rates,  ranging  from  100/s  to  1E9  Is.  The  lowest  strain  rate  produced  a  maximum  stress  of 
about  2.4GPa  while  the  highest  strain  rate  produces  a  maximum  stress  of  over  3.4GPa. 
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These  curves  were  developed  using  the  Johnson-Cook  constitutive  equation.  The 
critical  stress  is  the  maximum  stress  on  the  curve,  and  the  critical  equivalent  plastic  strain 
is  the  strain  that  occurs  at  the  maximum  stress  value.  To  find  the  critical  strain,  one  plots 
the  critical  plastic  strain  against  the  strain  rate  and  fits  performs  a  curve  fit.  This  was 
done  by  Meador,  and  his  resulting  plot  is  shown  below  in  Figure  21. 


Johnson-Cook  Plastic  Strain  at  Max  Flow  Stress  vs  Plastc  Strain  Rate 


Figure  21:  Critical  Plastic  Strain  Curve  Fit  [28] 

The  resulting  curve  fit  is  shown  in  Equation  2.30. 

2  30 

£CritO<  y,  0  =  APSe(x,  y,  t)Bps  +  CPS 

Where  ecrit  is  the  critical  strain,  e  is  the  strain  rate,  Aps,  Bps,  and  Cps  are  curve  fit 
constants. 
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Table  3  shows  the  result  of  the  curve  fit. 


Table  3:  Coefficients  of  Plastic  Strain  [4] 


Coefficient 

Value 

Aps 

2.247E-2  MPa 

Bps 

-5.516E-2 

Cps 

6.044E-3  MPa 

Note  that  unlike  the  Johnson  and  Cook  Failure  Criteria,  damaged  material  remains 


in  the  simulation. 

2.7.3  von  Mises  Stress  Criterion 

With  the  von  Mises  stress  failure  criterion,  there  is  a  maximum  stress  value  that 
that  once  reached,  counts  the  material  as  failed.  This  criterion  was  developed  by  Huber 
using  work  done  by  Hale.  Hale  found  in  his  research  that  there  were  dominant  strain  rates 
for  certain  velocity  ranges.  [3]  A  maximum  stress  value  from  Johnson  and  Cook’s  model 
for  VascoMax300  was  determined  by  Hale  based  on  the  dominant  strain  rates.  Hale’s 
results  are  seen  in  Table  3  below: 


Table  4:  Maximus  Stress  at  Dominant  Strain  for  VascoMax  300  [3] 


Velocity  Range  (m/s) 

Dominant  Strain 

Rate 

Maximum  Stress 
(MPa) 

10-200 

1  X  105 

2,900 

300-622 

1  x  106 

3,000 

750-1530 

1  x  107 

3,130 

Based  on  Hale’s  results,  a  maximum  stress  value  of  3  GPa  was  set  and  compared 


to  the  von  Mises  stress  calculated  within  CTH.  As  with  the  strain  at  max  stress  criterion, 
while  the  material  is  marked  as  failed,  it  is  not  removed  from  the  simulation  as  it  is  with 
the  Johnson  and  Cook  failure  criterion. 
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2.8  Hydrocodes 

Due  to  the  complexity  and  the  sheer  magnitude  of  the  possible  shock  problems,  it 
is  of  great  help  to  have  numerical  solvers  that  can  simultaneously  solve  these  equations. 

A  hydrocode  is  an  extension  of  computational  fluid  dynamics  programs.  They  are 
particularly  geared  for  simulating  dynamic  events  involving  shocks.  There  are  two 
reference  frames  that  describe  the  way  materials  deform:  Eulerian,  or  Lagrangian.  [36] 

Eulerian  and  Lagrangian  systems  vary  in  how  they  spatially  represent  materials  in 
a  system.  Figure  22  gives  a  visual  representation  to  show  the  difference. 


Figure  22:  Lagrangian  and  Eulerian  Coordinate  Systems  [41] 

Consider  the  case  above  where  a  rectangle  is  accelerating  into  a  wall  and 

deforming.  In  the  Lagrangian  representation,  the  grid  points  are  attached  to  the  material 
and  move  with  the  material  as  the  simulation  progresses.  Note  how  meshing  is  only 
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needed  around  the  materials.  Because  the  cells  are  allowed  to  deform,  the  volume 
changes,  but  never  the  mass,  so  thus  the  density  is  variable  depending  on  if  the  cell  is 
compressed  or  stretched. 

In  the  Eulerian  case,  material  is  placed  inside  the  mesh  -  which  covers  an  entire 
area  of  action  -  and  the  grid  points  never  move.  Calculations  are  made  based  on  the  flow 
rates  of  materials  through  the  cells.  Thus,  the  volume  of  the  cell  never  changes,  but  the 
amount  of  mass  present  in  the  cell  can,  and  thus  also  the  density.  [36] 

The  equations  used  for  each  approach  differ  slightly  due  to  the  differing  states. 
Below  are  the  Lagrangian  and  Eulerian  conservation  equations  used  by  hydrocodes.  [36] 
Lagrangian  equations  for  conservation  of  mass,  momentum,  and  energy 
respectively: 


Dp  dvt 

M+pa^-° 

2.31 

Dv;  1  do  a 

=  f  + 

Dt  1  p  dxj 

2.32 

De  1  d  ,  . 

—  =  fjVj  +  — - — (OiiVi) 

Dt  n  1  pdXjK  lJ  lJ 

2.33 

Eulerian  equations  for  conservation  of  mass,  momentum,  and  energy  respectively: 

dp  d 

Tt+wyd=o 

2.34 

dvt  dvt  _  1  doji 

dt  Vi  dxj  1  p  dxj 

2.35 

de  de  Id,. 

-+Vi—=fiVl+-—(ailVi) 

2.36 
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Where: 

Vi  is  the  velocity 

p  is  the  density 

Oji  is  the  stress  tensor 

e  is  the  specific  total  energy 

fi  =  the  body  force 

xt  is  the  Eulerian  coordinate  axis 

i,  and  j  are  summation  indices 

t  is  time 


It  is  further  noted  that  the  total  energy  e  can  be  split  into  the  kinetic  energy  and 
the  internal  energy  as  seen  in  the  equation  below: 


e  =  -v.-Vi  +  E 
2 


2.37 


Where  E  is  the  specific  internal  energy. 

With  this  addition,  the  conservation  of  energy  equation  can  be  written  for  the 
Lagrangian  and  Eulerian  cases  respectively  as: 


DE  P  Dp  1 

— —  —  — —— - 1 - SjjEjj 

Dt  p2  Dt  p 


2.38 

2.39 


dE  dE  P  (dp  dp\  1 

w+v‘d^=p{w+v‘d 

Where  is  the  stress  deviators,  £i;-  is  the  strain  rates,  and  P  is  the  hydrostatic 


pressure. 

The  Lagrangian  and  Eulerian  equations  are  interchangeable  by  use  of  the 
substantial  derivative,  where: 

D  d  d 

—  =  —  +  vt  - —  2.40 

Dt  dt  dxt 
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The  program  that  is  used  in  this  research,  CTH,  is  based  primarily  in  the  Eulerian 
system,  but  it  does  use  the  Lagrangian  conservation  equations  before  mapping  them  back 
over  to  the  Eulerian  system.  This  will  be  further  addressed  in  the  following  sections. 

Both  sets  of  equations  calculate  solutions  at  grid  points  of  meshes  that  are  defined 
in  the  code. 

Both  systems  have  different  advantages  and  disadvantages.  The  Lagrangian 
system  can  show  the  shape  of  the  deformed  materials  much  more  accurately,  and  can 
work  more  easily  with  constitutive  equations  and  failure  criteria,  because  it  is  able  to 
record  the  time  history  of  the  material.  However,  if  events  are  too  dynamic  and  the  cells 
become  extremely  distorted,  the  Lagrangian  system  can  struggle  with  making  said 
calculations.  The  Eulerian  system  is  much  better  in  this  case  and  is  able  to  allow  for  these 
extreme  cases  to  occur.  However,  the  material  position  is  only  vaguely  known,  because 
once  a  material  enters  a  cell,  it  is  considered  present  in  the  entire  cell.  Thus  the 
calculations  for  the  constitutive  equations  and  failure  criteria  are  doable,  but  not  as 
accurate  in  a  Lagrangian  system.  Additionally,  since  materials  are  prone  to  mixing,  the 
boundary  conditions  can  be  a  bit  trickier  in  an  Eulerian  system. 

2.9  CTH 

The  specific  hydrocode  used  in  this  research  attempts  to  bring  together  the  best  of 
both  worlds.  It  was  developed  by  Sandia  National  Laboratories  to  model  ballistic 
problems,  so  excels  at  modeling  strong  shock  wave  physics.  It  is  an  Eulerian-based 
system,  but  uses  a  Lagrangian  step  in  its  solving.  [39]  The  solving  process  uses  explicit 
integration  based  on  a  time  step  In  the  process  of  making  calculations,  first  CTH  takes  a 
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Lagrangian  step  and  solves  the  Lagrangian  conservation  equations  for  the  given  time 
step,  then  the  deformed  mesh  is  re-mapped  back  to  the  Eulerian  mesh.  A  new  time  step  is 
then  calculated  to  satisfy  the  stability  of  the  solution.  [25] 

2.9.1  Lagrangian  Step 

During  the  Lagrangian  step,  the  mesh  is  distorted  to  follow  the  motion  of  the 
material  in  Lagrangian  using  the  Lagrangian  conservation  equations.  The  conservation  of 
mass  is  already  fulfilled  since  mass  cannot  exit  the  cell.  Energy  and  momentum  are 
calculated  using  finite  volume  representations.  Care  must  be  taken  to  select  a  time  step 
that  is  small  enough  to  prevent  a  wave  from  completely  crossing  a  cell  in  one  step.  Mass, 
volume,  momentum,  and  energy  are  all  calculated.  In  the  event  that  the  conservation 
equations  cannot  be  satisfied,  CTH  will  always  choose  to  satisfy  the  momentum  equation. 

2.9.2  Remap  Step 

Calculations  made  in  the  Lagrangian  step  are  transferred  and  new  calculations  are 
made  to  align  them  with  the  Eulerian  mesh.  Change  in  volume  from  the  old  cell  to  the 
new  cell  and  where  the  changing  masses  end  up  are  determined  by  an  interface  tracking 
algorithm  set  in  CTH.  Mass  and  internal  energy  is  shifted  between  old  and  new  cells. 
Using  data  from  the  tracking  algorithm,  momentum  and  kinetic  energy  is  divided 
between  cells  based  on  the  movement  of  mass.  After  the  re-map  step  is  complete,  the 
next  time  step  is  calculated. 

2.9.3  Boundary  Conditions 

Boundary  conditions  are  also  an  essential  part  of  bringing  a  numerical  model 
together.  CTH  has  6  different  types  of  acoustic  wave  oriented  boundary  conditions  and 
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two  of  them  are  used  in  this  model.  The  mesh  that  is  set  up  is  surrounded  by  a  secondary 
outer  mesh  consisting  of  “ghost  cells.”  [39]  These  cells  are  used  purely  for  numerical 
calculations,  and  the  set  boundary  conditions  determine  how  mass,  momentum,  and 
energy  fluxes  across  boundaries. 

The  Type  0  boundary  is  the  symmetry  boundary  condition,  which  models  rigid 
boundaries  and  symmetrical  geometry.  This  is  particularly  handy  for  cylindrical  and 
spherical  geometries.  Velocity  between  the  ghost  cell  and  the  interior  mesh  cell  is  set  to  0 
and  the  kinetic  energy  is  converted  to  internal  energy;  there  is  no  flux  across  the 
boundary. 

The  Type  1  boundary  is  the  sound  speed  based  absorbing  boundary  condition, 
which  simulates  and  infinite/semi-infinite  medium.  Mass  is  allowed  to  flow  in  and  out  of 
the  mesh.  This  type  is  used  on  the  bottoms  of  the  x  and  y  meshes  in  our  scenario,  where 
the  material  is  situated. 

The  Type  2  boundary  is  the  outflow  boundary  condition.  It  allows  mass  to  leave 
the  system.  When  mass  flows  into  the  ghost  cells,  the  ghost  cells  are  automatically  set  to 
void.  This  is  the  condition  places  on  the  top  of  the  x  and  y  axis,  where  there  is  no  material 
initially. 

The  Type  3  boundary  condition  is  the  inflow  boundary  condition,  where  mass  can 
enter  the  interior  of  the  mesh  from  the  outside.  The  properties  of  the  incoming  mass  are 
extrapolated  from  adjacent  interior  cells. 

The  Type  4  Boundary  condition  is  the  outflow  extrapolation  boundary  condition. 

It  linearly  extrapolates  the  boundary  pressure  from  the  interior  of  them  mesh.  If  the 
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extrapolated  pressure  turns  out  to  be  positive,  then  it  is  set  to  zero.  Mass  can  flow  in  and 
out  of  the  mesh.  This  condition  is  good  for  letting  parts  of  the  problem  flow  out  while 
maintaining  pressure, 

The  Type  5  boundary  condition  is  similar  to  the  Type  0  condition.  It  is  symmetric, 
but  mass  is  not  allowed  to  flow  in. 

Type  6  boundary  is  periodic.  This  one  is  particularly  unique  because  it  has 
specific  processor  requirements  (at  least  2  processors  must  be  in  the  periodic  direction), 
and  the  boundaries  must  be  set  to  be  symmetrical,  so  the  top  and  the  bottom  must  both  be 
periodic.  It  is  used  for  modeling  infinite  or  semi-infinite  materials  with  an  explicit 
micro  structure. 

While  there  are  many  options  for  boundary  conditions  in  CTH,  only  type  1  and 
type  2  are  used  in  our  simulation  to  simulate  the  infinite  air  medium  and  to  allow  the 
slipper  to  exit  the  mesh. 

2.9.4  Mixed  Materials 

What  does  CTH  do  when  more  than  one  material  occupies  a  cell?  CTH  has  a  few 
different  options  for  how  to  deal  with  it.  With  mixed  cell  option  1  (Mix  1),  the  strength  in 
the  cell  is  a  weighted  average  based  on  the  volume  fractions  of  the  materials  in  the  cells. 
Mix  5  sets  the  yield  strength  to  0  when  materials  mix.  Mix  3  is  like  mix  1  in  that  it  bases 
the  yield  strength  on  the  mix  of  the  materials,  but  divided  by  the  sum  of  the  volume 
fraction  of  the  materials  that  support  shears.  This  is  what  generally  gives  the  best  results 
and  what  is  used  in  our  model.  [39]  If  the  reader  is  curious  what  happened  to  mixes  2  and 
4,  they  were  options  in  previous  generations  of  CTH  but  have  since  been  removed. 
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2.9.4  Material  Interface 

There  are  two  options  for  how  material  will  behave  when  it  comes  into  contact 
with  another  material.  The  options  are  to  allow  for  a  sliding  interface  (“slide”),  where 
there  is  no  friction,  or  to  call  for  a  “no  slide”  condition,  where  the  materials  are 
essentially  “welded  together.”  [39]  The  materials  will  not  move  forward  until  the 
averaged  fracture  pressure  of  the  materials  in  the  cells  has  been  reached.  The  reason  for 
this  very  on/off  approach  to  materials  interacting  is  that  in  strong  shock  wave  scenarios 
the  normal  stresses  are  usually  much  stronger  than  shear  stresses,  so  typically  fine  tuning 
them  is  rarely  necessary.  Sadly,  this  means  that  the  friction  properties  for  sliding  are  non- 
negotiable.  However,  using  the  no  slide  condition  gives  an  upper  bound  for  wear  because 
far  more  damage  to  the  slipper  must  take  place  for  the  slipper  to  even  move,  and  the 
slipper  is  continually  taking  additional  damage  as  it  is  sliding.  Dave  Huber  did  work 
comparing  the  slide  and  no  slide  conditions,  and  as  expected,  wear  rates  using  the  no 
slide  condition  were  always  much  higher  than  those  using  the  slide  condition.  [29] 

2.9.5  Tracers 

Tracer  points  can  be  set  up  within  CTH  to  track  the  progress  of  variables. 
They  are  set  at  initial  Lagrangian  points  and  record  requested  data,  such  as  the  velocity, 
the  pressure,  the  yield  stress,  the  temperature,  and  any  other  data  the  user  may  want.  The 
user  can  even  define  variables  and  have  CTH  perform  calculations  during  the  simulation 
and  output  them  via  tracer  points.  The  data  collected  is  what  allows  for  post-processing 
and  calculation  of  wear  values  to  occur. 


47 


To  track  what  exactly  is  happening  as  CTH  progresses,  let  us  examine  the  process 
of  the  failure  of  a  cell  using  the  Johnson  and  Cook  failure  criterion.  At  the  initial  stage,  a 
velocity  is  imposed  on  the  material  along  with  a  given  time  step.  The  first  Lagrangian 
step  occurs  and  the  mesh  deforms  while  state  variables  and  flow  stress  are  calculated, 
along  with  the  plastic  strain  and  the  plastic  strain  at  failure  for  the  damage  criteria.  The 
mesh  is  then  re-mapped  back  to  its  original  position,  and  the  calculated  variables  are  also 
shifted  into  the  Eulerian  system.  The  failure  criteria  is  checked  -  has  it  reached  its  critical 
value?  If  not,  the  next  time  step  is  calculated  and  the  process  repeats  itself  until.  When 
using  the  Johnson  and  Cook  failure  criterion,  if  the  critical  value  has  been  reached,  the 
strength  of  the  material  in  that  cell  is  set  to  0,  so  it  can  no  longer  support  other  material 
and  behaves  essentially  like  a  fluid.  If  not,  then  the  cell  will  fail  based  on  when  a  fracture 
pressure  of  the  materials  is  reached. 

2.10  Defining  Wear 

Our  research  is  interested  in  bringing  together  all  of  these  pieces  to  determine 
wear.  Wear  here  is  the  volume  of  material  removed  from  the  slipper.  Wear  rates  are  a 
ratio  of  the  volume  of  material  removed  to  the  distance  slid.  Additionally  wear  rates  per 
unit  area,  are  used,  which  is  the  wear  area  divided  by  the  distance  slid. 

Wear  is  defined  at  the  volume  of  material  removed.  To  obtain  values  for  the  wear 
rate  in  mm7mm,  one  runs  a  simulation  and  determines  the  number  of  failed  cells  based 
on  whatever  fracture  criterion  is  being  used,  uses  that  number  to  calculate  a  volume  and 
divides  it  by  the  distance  slid. 
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2.11  Summary  of  Theory 

In  this  section,  the  background  behind  the  equations  being  used  in  this  analysis 
has  been  reviewed,  including  wave  mechanics  through  elastic,  plastic,  and  shock  waves, 
fundamental  conservation  equation  reviewed,  Lagrangian  and  Eulerian  states  explained, 
equations  of  states,  constitutive  equations,  and  failure  criterion  reviewed.  Following  this, 
a  brief  overview  of  how  CTH  applied  these  building  blocks  in  its  solving  procedure  was 
given.  Now  that  the  problem  background  and  theory  have  been  addressed,  let  us  move  on 
to  previous  wear  results  that  this  research  attempts  to  match. 
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III.  Comparative  Data 


3.1  Chapter  Overview 

This  research  is  focused  around  comparing  wear  results  to  experimental  and  other 
modeled  data.  Major  Buentello  has  been  working  on  an  ABAQUS  FEA  model  that 
models  the  same  scenario  as  the  CTH  model.  His  model  considers  a  few  different 
collision  scenarios,  and  is  able  to  incorporate  Coulomb’s  friction  model  (using  a 
coefficient  of  friction  determined  in  previous  research  by  Montgomery)  and  uses  Mie- 
Gruneisen  EOS  rather  than  the  SESAME  tables.  [10]  As  such,  he  also  has  additional 
thermal  considerations  that  are  outside  the  scope  of  the  CTH  model.  The  idea  in  this 
research  is  to  see  if  some  of  the  refinement  incorporated  into  the  finite  element  approach 
is  really  necessary.  Therefore,  the  approach  used  herein  is  to  make  use  of  a  simple  plane 
strain  model  with  corrections  and  evaluate  wear  volume  with  the  defined  damage 
relations.  The  end  result  will  then  be  compared  without  a  thermal  environment  to  the 
finite  element  solution. 

3.2  2008  Experimental  DADS  Data 

The  HHTSS  data  this  research  is  attempting  to  compare  to  come  from  a  2008 
experiment.  The  set-up  for  the  particular  experiment  can  be  seen  in  Figure  23. 
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Figure  23:  Rocket  Stages  for  the  2008  HHSTT  Experiment  [3] 

It  consisted  of  the  front  forebody  sled,  which  carried  the  experiment  payload,  and 

three  pusher  sleds,  which  carried  rockets  that  propelled  the  experiment  up  to  l,530m/s. 

The  entire  experiment  occurred  in  8.14  seconds  and  over  a  distance  of  5,816  m.  [38]  The 

front  two  sleds  were  fitted  with  slippers  made  from  VascoMax300,  but  the  back  two 

sleds’  slippers  were  made  from  AISI4130  Steel.  One  of  the  VascoMax300  slippers  from 

stage  3  was  recovered  and  it  is  the  one  focused  on  in  the  research.  A  diagram  of  the 

slipper  can  be  seen  in  Figure  24. 


Figure  24:  Slipper  and  Rail  dimensions  [3] 
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From  technical  drawings,  the  un-used  slipper  is  assumed  to  have  a  nominal 
thickness  of  14.7  mm.  The  area  of  the  slipper  subjected  to  wear  damage  was  determined 
to  be  103.5  mm  ,  and  the  width  of  the  slipper  was  determined  to  be  10.1  mm.  [38]  The 

-j 

total  volumetric  wear  was  determined  to  be  10.52  cm  ,  or  3.49%  of  the  initial  top  volume 
(300.9cm  ).  The  velocity  profile  of  the  slipper  is  shown  below  in  Figure  25. 


Sled  C.G.  Velocity 


Figure  25:  Horizontal  Sled  Velocity  of  a  stage  3  slipper  from  the  2008  HHSTT  Experiment  [3] 
Though  one  would  like  to  imagine  the  entire  trip  is  a  smooth  ride,  there  are  some 

complications  to  consider.  For  one  thing,  the  slipper  does  not  stay  in  consistent  contact 

with  the  rail.  Figure  26  shows  an  incriminating  shot  of  the  slipper  during  a  similar 

experiment  in  2003. 
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(a)  80X-G1  Test  Mission,  April  2003,  V  =  2,853  m/s. 


(b)  80X-G1  Test  Mission,  April  2003.  V  =  2,861  m/s. 


Figure  26:  2003  HHSTT  Experiment  showing  Slipper  Pitch 

As  the  picture  shows,  the  slipper  is  pitched  forward,  so  the  majority  of  the  wear  is 

occurring  at  the  front  of  the  slipper.  Additionally,  due  to  the  aerodynamics  of  the 
experiment,  lift  is  generated,  causing  the  slipper  to  bounce  on  the  rail. 
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Figure  27:  Discretized  Horizontal  Velocity,  Vertical  Position,  and  Vertical  Velocity  of  a  Stage  3 

Slipper  [3] 

This  figure  shows  the  vibration  while  the  slipper  runs  along  the  track.  This  means 
that  the  bottom  of  the  slipper  is  not  in  contact  with  the  rail  100%  of  the  time.  This  has 
already  affected  research  with  Meador  where  he  discovered  that  due  to  the  bouncing  the 
thermal  distribution  coefficient  needed  to  be  adjusted.  [28]  It  is  a  major  consideration  that 
Buentello  is  attempting  to  account  for  in  his  ABAQUS  model. 
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3.3  3D  Abaqus  Model 


Major  Buentello  has  created  a  model  in  ABAQUS  to  attempt  to  model  the  wear.  It 
models  the  same  scenario  as  the  CTH  model  -  a  slipper  colliding  with  single  asperity. 
His  model  uses  the  Johnsons  and  Cook  flow  stress  model,  implements  the  Mie- 
Griineisen  EOS  as  opposed  to  SESAME  tables,  and  implements  a  contact  algorithm 
between  the  two  materials  and  uses  Coulomb’s  Friction  Model  as  seen  below  in 
Equation  3.1. 

Ff<fiFn  3.1 

Where 

Ff  =  Friction  force 
p  =  Friction  coefficient 
Fn  =  Normal  force 

The  friction  coefficient  used  is  based  on  Montgomery’s  previous  work.  [10] 

n 

Vd  =  ^  WuAnACnDSn%Contn  3.2 

71  =  1 


Where: 

Vd=  Damage  volume 
WuAn=  Wear  per  unit  area 
Ac  =  Area  of  contact 
Ds  =  Sliding  distance 

%Cont  =  Percent  of  contact  time  the  slipper  had  with  the  rail. 
The  results  are  shown  below  in  Table  5. 
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Table  5:  ABAQUS  Volume  Damage  Calculations  [38] 


Avg. 

Vel 

(m/s) 

WuAi 

(mmA3/mmA2) 

width 

(mm) 

Avg 
Worn 
length 
li  (mm) 

Sliding 

Distance 

Di 

% 

Contact 

Predicted 

Vol. 

Damage 

Vi 

Measured 

Vi 

Stage 

1 

175 

1.67E-6 

101 

70 

323000 

.3063 

1168.1 

6716.5 

Stage 

2 

500 

1.81E-6 

101 

35 

1016000 

.3063 

1993.9 

1684.2 

Stage 

3 

1100 

2.65E-6 

101 

22.5 

1356000 

.3063 

2501.2 

1010 

Desa 

1 

1250 

2.78E-6 

101 

10 

1265000 

.3063 

1086.0 

252.5 

Desa 

2 

750 

2.2E-6 

101 

25 

1680000 

.3063 

2858.5 

909 

Totals 

5816000 

9607.8 

10572 

The  Wua  numbers  are  generated  from  Buentello’s  ABAQUS  model  using  a  2um 
asperity  for  the  collision.  The  width  is  known  from  technical  drawings  and  the  sliding 
distance  and  is  known  from  the  DADS  data.  The  percent  contact  is  assumed  to  be  0.3063 
from  previous  work.  The  predicted  volume  damage  was  calculated  using  equation  3.2. 
The  measured  volume  came  from  calculations  made  by  Buentello.  Looking  over  work 
done  by  Hale,  Buentello  noted  that  there  were  5  general  worn  areas  on  the  slipper  and 
was  able  to  identity  which  area  would  have  worn  away  based  on  each  acceleration  or 
deceleration  stage.  This  is  illustrated  in  Figure  28. 
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Initially  End  of  Stage  1  End  of  Stage  2  End  of  Stage  3  Deceleration 


Figure  28:  Wear  Areas  based  on  Acceleration  State  [38] 

From  this  observation,  he  calculated  the  observed  removed  volume  for  each  stage. 

The  total  wear  calculated  from  this  method  closely  matched  with  the  total  wear  removed 

from  the  slipper  (10.52  cm  compared  to  10.57  cm  ).  This  research  attempts  to  compare 

its  results  to  the  Stage  2,  500mps  wear  area  rate  and  wear  volume. 


3.4  Differences  between  ABAQUS  and  CTH 

Previously  the  workings  of  CTH  were  described,  but  how  exactly  are  the 
differences  from  the  ABAQUS  code  used?  The  most  obvious  difference  is  that  CTH  is  an 
Eulerian-based  system,  but  ABAQUS  is  Lagrangian  based  system.  This  means  CTH  has 
an  advantage  at  higher  velocities  because  its  mesh  won’t  distort  too  much  to  be  able  to 
perform  calculations,  but  ABAQUS  had  better  defined  boundary  conditions  between 
materials.  While  CTH  uses  acoustic  wave  boundary  conditions,  ABAQUS  boundary 
conditions  are  based  on  stress  and  displacement.  Additionally  they’re  more  easily  defined 
since  its  mesh  is  attached  to  materials  separately,  and  materials  aren’t  sharing  as  in  CTH. 
Not  only  will  it  more  accurately  predict  specific  shapes  of  material  deformation,  but  it 
can  also  have  a  defined  friction  force  between  materials. 
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Figure  29:  ABAQUS  Model  showing  damage  occurring  to  the  slipper  [38] 

Thus,  ABAQUS  is  better  able  to  take  into  account  adhesion  and  abrasion  wear, 

while  CTH’s  ability  to  add  it  in  is  crude  at  best.  Additionally  the  way  the  two  programs 

make  then  calculations  differ.  CTH  uses  finite  volume  calculations  by  explicitly 

integrating  at  the  edges  of  the  cells  to  get  its  data.  ABAQUS  uses  finite  difference 

methods  to  get  its  data. 

3.5  Wolfson  Experimental  Data 

Wolfson  in  his  experiments  was  able  vary  the  sliding  materials  he  used,  velocities, 
normal  forces,  track  conditions,  and  bearing  pressures.  As  such  not  all  of  his  experiments 
are  applicable  to  this  study,  but  of  the  60  runs,  2  of  them  are  comparable.  They  use 
stainless  steel  sliders  on  a  bare  steel  track.  The  wear  was  originally  in  units  of  in./ft,  but 
by  multiplying  the  results  by  the  slider  area  and  converting  units  the  data  is  comparable 
to  this  study’s  results.  The  conversion  is  shown  in  Equation  3.3. 

Wearrate  ^woifson^n  3.3 

Where  Wwoifson  is  the  measured  wear  rate  and  Anis  1  in.2  (645.16mm2) 

The  values  used  are  seen  in  Table  6. 
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Table  6:  Wolfson  Data  Conversion 


Sliding  Velocity 

Average  Wolfson  Wear 
Rate 

Converted  Average 
Wear  Rate 

ft/s 

m/s 

in./ft 

mm/mm 

mnrVmm 

825 

251.46 

9.50E-6 

7.92E-7 

5.11E-4 

2500 

762 

7.50E-6 

6.25E-7 

4.03E-4 

3.6  Previous  AFIT  Work 

Both  Chad  Hale  and  Dave  Huber  have  done  previous  work  that  this  research  will 
compare  against.  Chad  Hale  developed  a  global  and  a  local  submodel  finite  element 
analysis.  The  results  from  his  local  submodel  are  used  to  compare  to  our  results.  He  ran 
12  different  target  velocities  using  a  6um  asperity.  His  model  dimensions  can  be  seen  in 
Figure  30. 


Sliding  Direction 


6  pm  Asperity 

\ 


-  75  pm 

Slipper  Submodel 


75 


pm 


Deformable  Rail 


Figure  30:  Hale's  Plain  Strain  ABAQUS  Model  Dimensions 
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Hale  uses  the  Mie-Gruneisen  EOS,  the  Johnson  and  Cook  constitutive 


Viscoplastic  model,  and  also  uses  the  Johnson  and  Cook  failure  criterion.  The  results 
from  his  work  are  shown  below  in  Table  7. 


Table  7:  Hale’s  Wear  Volume  Rate  Results 


Velocity 

Hale  (FEA) 

(m/s) 

mmA3/mm 

10 

2.63E-05 

20 

3.44E-05 

40 

3.79E-05 

100 

4.88E-05 

200 

8.42E-05 

300 

1.20E-04 

500 

2.00E-04 

622 

2.13E-04 

750 

1.93E-04 

1000 

1.87E-04 

1250 

1.67E-04 

1530 

1.37E-04 

Huber  used  a  model  near  identical  to  the  model  used  in  this  research.  The 
difference  between  his  work  and  this  work  is  the  use  of  the  Johnson  and  Cook  Failure 
Criterion,  and  the  ability  to  pinpoint  maximum  wear.  Huber’s  work,  simulations  were  run 
using  the  following  equation: 


tsim 


(1-1  )ra 

V  slide 


3.4 


Where  tsim  is  the  simulation  time,  ra  is  the  asperity  radius,  and  vsiide  is  the  sliding 
velocity. 

This  equation  was  used  so  that  the  slipper  would  run  over  the  maximum  height  of 
the  asperity,  however  it  was  found  in  this  research  that  the  time  for  maximum  wear  varied 
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based  on  failure  criteria  and  slipper  size.  The  data  used  from  Huber’s  work  is  shown  in 
Table  8. 


Table  8:  Huber's  Wear  Volume  Rate  Results 


Velocity 

Slide  Line  Strain  at  omax 

Slide  Line  von  Mises 

stress 

m/s 

mmA3/mm 

mmA3/mm 

200 

7.49E-05 

3.71E-05 

300 

8.34E-05 

6.16E-05 

400 

1.05E-04 

1.09E-04 

500 

1.45E-04 

1.70E-04 

600 

1.70E-04 

2.04E-04 

700 

1.86E-04 

2.42E-04 

800 

2.06E-04 

2.79E-04 

900 

2.20E-04 

3.03E-04 

1000 

2.32E-04 

3.21E-04 

1100 

2.44E-04 

3.33E-04 

1200 

2.58E-04 

3.44E-04 

1300 

2.62E-04 

3.63E-04 

1400 

2.64E-04 

3.61E-04 

1500 

2.67E-04 

3.45E-04 
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IV.  Analysis  and  Results 


4.1  Goal 

The  aim  of  this  research  is  to  implement  a  2D  plane  strain  model  in  CTH  using 
the  Johnson  and  Cook  fracture  criterion  to  represent  a  3D  wear  event,  to  validate  it  by 
comparing  it  to  previous  studies’  wear  results. 


4.2  Plane  Strain  Model 

It  is  of  interest  to  use  a  2D  plane  strain  model  because  of  its  simplicity. 
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Figure  31:  Plane  Strain  Diagram  [X] 

By  placing  the  model  in  a  2D  coordinate  system,  the  model  completely  ignores 
strains  in  the  z  direction  and  shear  strains  in  the  xz  and  yx  planes. 

The  2D  plane  strain  model  consists  of  a  0.0125cm  by  0.0649cm  slipper  sliding 
along  a  0.0200cm  by  0.850cm  rail  with  a  constant  velocity.  The  total  model  is  seen  in 
Figure  32. 
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Figure  32:  CTH  Plane  Strain  Model 

The  slipper  collides  with  an  asperity  placed  in  front  of  it,  which  range  in  size  from 
6um  to  2um.l295  cells  are  traced  throughout  the  simulation,  and  recorded  are  the 
position,  yield  stress,  cell  volume  fraction  (fraction  of  material  volumes  present  in  each 
cell),  strain  rate,  and  optionally  the  Johnson  and  Cook  damage  coefficient.  In  this 
research  velocities  varied  between  200m/s  and  1500m/s. 

To  get  an  idea  of  how  the  material  collision  looks  in  CTH,  Figures  33-36  show 
the  progress  of  the  material  flow  over  time. 
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Figure  33:  Material  Flow  Stage  1 


685  690  695  700  705  710  715 
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Figure  34:  Material  Flow  Stage  2 
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Figure  35:  Material  Flow  Stage  3 

Materials  at  1 ,41e-08  seconds 
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Figure  36:  Material  Flow  Stage  4 


It  is  seen  that  the  boundaries  between  the  materials  aren’t  exactly  crisp,  as  is 
characteristic  of  an  Eulerian  system.  Deformation  is  not  really  occuring,  but  the  materials 
are  flowing  together. 

As  the  impact  occurs,  the  pressure  in  the  materials  increases  dramatically.  Figures 
37-40  show  how  this  progresses  in  the  simulation. 
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Figure  37:  Pressure  Development  Stage  1 

Pressure  at  7.00e-09  seconds 
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Figure  38:  Pressure  Development  Stage  2 
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Pressure  at  1 .06e-08  seconds 


660  680  700  720  740 


X  (pm) 

Figure  39:  Pressure  Development  Stage  3 


Pressure  at  1 ,47e-08  seconds 


660  680  700  720  740 
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Figure  40:  Pressure  Development  Stage  4 
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As  the  figures  show,  pressure  propagates  progressively  in  a  wave  pattern  in  both 
the  slipper  and  the  asperity,  with  the  less  intense  elastic  waves  traveling  more  quickly 
than  the  stronger  plastic  waves.  As  time  continues  the  strength  begins  to  dissipate, 

The  materials  for  the  slipper  and  rail  are  Vascomax300  and  1080  steel 
respectively,  to  mirror  what  is  actually  used  at  HHSTT.  The  properties  of  these  materials 
are  listed  below. 


Table  9:  VascoMax300  and  1080  Steel  Material  Properties 


Property 

VascoMax  300 

1080  Steel 

Melt  Temperature,  Tmen 

1,685K 

1,670  K 

Density,  p 

8,  000  kg/m3 

7,800  kg/m3 

Poisson’s  Ratio,  v 

0.283 

0.270 

Modulus  of  Elasticity,  E 

180.7  GPa 

202.8  GPa 

Hardness,  H 

200  GPa 

- 

The  Johnson  and  Cook  plasticity  model  and  failure  criterion  are  also  used  in  the  model, 
with  constants  as  shown  below. 

a  =  (A  + Ben)(l  + C\ns*)(l-T*m)  2.27 


Table  10:  Johnson  and  Cooke  Constitutive  Model  Material  Constants 


Constant 

VascoMax  300 

Steel  1080 

A 

2.1E10 

0.7E10 

B 

0.124E10 

3.6E10 

C 

0.03 

0.17 

M 

0.8 

0.25 

N 

.3737 

0.6 

D  = 


de 

£pf  ( T,p,£ ,  £) 


2.28 


£pf  =  [Dl  +  D2exp(-D3^)]  [D4  +  ln(max(l,  £))][!  +  DST*] 


2.29 
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Table  11:  Johnson  and  Cook  Failure  Criterion  Material  Constants 


Constant 

VascoMax  300 

D1 

-0.09 

D2 

0.27 

D3 

-0.48 

D4 

0.014 

D5 

3.870 

One  thing  to  note  about  the  Johnson  and  Cook  Failure  Criterion  constants  is  that 
D3  is  negative  in  CTH.  The  way  the  equation  in  arranged  in  ABAQUS,  D3  is  positive. 

4.3  Wear  Calculation  Road  Map 

Wear  is  defined  at  the  volume  of  material  removed.  To  obtain  values  for  the  wear 

-1 

rate  in  mm7mm,  one  runs  a  simulation  and  determines  the  number  of  failed  cells  based 
on  whatever  fracture  criteria  being  used,  uses  that  number  to  calculate  an  area  and  divides 
it  by  the  distance  slid.  This  number  is  then  multiplied  by  a  spherical  correction  factor 
(SCF).  Specifics  for  how  this  is  done  follow. 

To  give  an  overview  of  what  is  to  be  presented,  first,  how  the  damaged  areas  for 
each  failure  criteria  will  be  presented  and  then  used  to  calculate  wear  per  unit  width 
(mm7mm)  in  section  4.4.  Results  for  these  calculations  are  presented  graphically.  The 
next  topic  introduced  is  how  to  calculate  the  SCF  in  section  4.5,  which  entails  how  to 
calculate  single  asperity  wear.  Results  for  the  single  asperity  wear  and  the  SCFs  are 
given.  Finally  in  section  4.6,  wear  volume  rates  (mm7mm)  are  given  and  results  are 
compared  to  prior  work. 
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4.4  Wear  Per  unit  Width  Calculaitions 


4.4.1  Plastic  Strain  at  Max  Stress  Criterion 

In  CTH,  the  strains  and  strain  rate  data  are  extracted  from  the  simulation.  Then 
the  critical  strain  equation  is  applied  to  the  strain  rate  data,  and  if  the  strain  is  over  the 
critical  strain  value,  then  the  cell  is  marked  as  failed  with  a  “1,”  non- failed  cells  are  set  to 
“0.”  The  table  of  failed-non-failed  cells  is  then  multiplied  by  the  volume  fraction  of  the 
slipper,  so  that  only  the  slipper  failed  cells  are  counted  and  not  any  rail  failures.  The 
finally  result  is  summed  up  to  then  apply  to  an  area. 

4.4.2  von  Mises  Stress  Criterion 

Based  on  Hale’s  results,  a  maximum  stress  value  of  3  GPa  was  set  and  compared 
to  the  von  Mises  stress  calculated  within  CTH.  [29]  Once  again,  cells  that  reached  the 
critical  value  are  set  to  1  and  all  other  cells  are  set  to  0.  The  failed  slipper  cells  are 
isolated  by  multiplying  by  the  slipper  volume  fraction.  The  final  result  is  summed  up  into 
a  single  value 

4.4.3  Johnson  and  Cook  Fracture  Criterion 

Johnson  and  Cook’s  Fracture  Criterion  was  addressed  previously,  and  is 
conveniently  calculated  within  CTH.  It  is  easier  to  use  as  a  criterion  because  it  is  already 
in  a  0-1  range,  allows  partially  failed  cells  to  be  added  into  the  calculation,  and  is  already 
only  tracking  the  slipper  material,  so  does  not  need  to  be  isolated  out  from  the  rail 
material  with  the  volume  fraction.  The  final  damage  numbers  in  the  simulation  are  simply 
summed. 
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4.3.4  Damage  Area 

Once  the  damage  sums  are  calculated  based  on  each  stress  criteria,  these  numbers 
are  converted  to  an  area  (mm")  by  multiplying  the  by  the  mesh  area  of  the  traced  cells.  A 
wear  rate  can  be  obtained  by  dividing  the  wear  area  by  the  distance  slid  to  have  a  final 
result  in  terms  of  mm7mm.  However,  to  make  comparison  with  3D  models,  the  wear  rate 
needs  to  be  in  mm7mm  units.  Thus  a  correction  factor  is  needed  to  apply  to  the  results. 

4.3.5  Putting  it  all  together 

To  obtain  the  wear  per  unit  width,  the  damage  area  is  calculated  depending  on 
what  failure  criteria  is  being  used,  and  then  divided  by  the  distance  slid.  [3] 

Ad 

W  -  —  4  1 

uw  Ds 

Where  Wuw  is  the  wear  per  unit  width,  Ad  is  the  damage  area  which  varies  depending  on 
the  failure  criterion  in  use,  and  Ds  is  the  sliding  distance  of  the  slipper. 

Previously  the  wear  per  unit  width  values  were  calculated  using  equation  3.4,  but 
by  plotting  the  wear  rate  against  time,  it  has  been  observed  that  this  time  is  not  when  the 
maximum  damage  to  the  slipper  has  occurred.  Work  has  been  done  in  ABAQUS  to  find 
the  time  when  the  maximum  damage  has  occurred  and  no  additional  damage  is  done  to 
the  slipper.  Similar  efforts  were  made  using  CTH,  but  the  results  are  inconclusive.  The 
specific  CTH  and  MATLAB  codes  used  can  be  found  in  Appendix  A,  B.l  and  B.2. 
Figures  41-43  show  some  examples  of  the  damage  area  per  time  results. 
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Figure  41:  6um  lOOOmps 
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Figure  42:  4um  1500mps 
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Figure  43:  lum:  500mps 

These  plots  show  different  asperity  sizes  and  different  collision  speeds.  No  two 
plots  are  exactly  alike.  When  the  maximum  damage  has  been  reached,  the  wear  area 
should  level  off  to  a  flat  line.  In  Figure  41,  the  von  Mises  damage  area  is  very  clearly 
exhibiting  this  behavior.  However,  the  strain  at  max  stress  and  Johnson  and  Cook  damage 
areas  do  not  level  off  and  are  continuing  to  increase.  Thus  it  is  thought  that  simulations 
will  need  to  be  run  for  longer  simulation  times  to  try  to  capture  maximum  damage  area 
for  the  other  two  criteria.  The  results  for  the  simulations  that  were  run  can  be  found  in 
Appendix  C. 

4.5  Spherical  Correction  Factor 

There  is  a  lot  that  can  be  done  with  the  wear  per  unit  width  data.  Chad  Hale  and 
Dave  Huber  were  able  to  compare  their  work  to  Wolfson’s  experimental  data.  But  how 
were  they  able  to  convert  their  data  from  mm'/mm  to  mm  ?  They  needed  a  way  to 
transition  their  results  from  a  2D  scenario  to  a  3D  scenario.  Hale  was  able  to  do  this  by 
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developing  a  conversion  factor,  which  in  this  research  is  referred  to  as  the  spherical 
correction  factor. 


6  pm  Radius 


(a)  Hemispherical  Asperity 

Hale  found  this  constant  by  finding  the  wear  of  a  single  asperity,  and  dividing  it 
by  the  wear  rate  per  unit  width  of  the  asperity  of  interest  (in  his  case,  6um). 


4.2 


Where  Wsa  is  the  single  asperity  wear  and  Wuw  is  the  wear  per  unit  width. 

This  was  done  for  multiple  heights  along  the  total  asperity  size  of  interest  for 
multiple  collision  velocities.  The  results  were  then  averaged  together  for  a  final  result. 

Hale  developed  this  constant  in  ABAQUS  and  Huber  used  it  in  his  work  in  CTH. 
However  since  the  number  can  have  a  very  large  impact  on  the  end  wear  results,  it  was  of 
interest  to  see  how  the  correction  factor  would  turn  out  if  it  were  calculated  in  CTH  using 
different  failure  criteria. 
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4.5.1  Single  Asperity  Wear  Rates 

To  be  able  to  get  the  correction  factor,  the  wear  rate  for  a  single  asperity  is 
needed.  The  single  asperity  wear  represents  the  wear  occurring  as  a  function  of  height  of 
the  asperity. 

asperity 

f  radius  ,  „  ,  a  o 

Wsa  =  2\  Wuw  (y)  dy  4.3 

Jo 

To  make  this  calculation,  wear  rates  per  unit  widths  were  calculated  for  different 
heights  along  the  asperity  -  for  a  6um  correction  factor,  heights  2,  4,  and  6  were  used. 
For  a  2um  asperity  factor,  0.5,  1,  and  2  microns  were  used.  All  of  these  wear  rates  were 
calculated  at  200,  500,  1000,  and  1500  m/s  in  contrast  to  Chad  Hale,  who  used  10,  100, 
500,  1000,  and  1530  m/s.  Lower  velocities  were  excluded  from  CTH  calculation  because 
CTH  does  not  handle  low  velocities  well  due  to  its  lack  of  friction  modeling.  These 
results  were  then  plotted  according  to  their  wear  rate  and  the  location  along  the  asperity 
laterally,  which  was  calculated  in  the  yz  plane  using  the  equation  below: 


Where  z  is  the  asperity  distance  in  the  z  plane,  r  is  the  asperity  radius,  and  y  is  the 
asperity  height  in  the  y  plane. 

This  was  done  and  the  results  are  as  follows  for  a  2um  asperity  using  the  strain  at 
max  stress  and  the  Johnson  and  Cook  Failure  criteria. 
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Table  12:  2um  Wear  per  Unit  Width  Data  for  Single  Asperity  Calculations  using  the  Strain  at  Max 

Stress  Failure  Criterion 


2um  Strain  at  Max  Stress  Failure 
Criteria 

Wear  rate 

oer  unit  widt 

h  (Wuw) 

200 

500 

1000 

1500 

Distance 

2  microns 

0.002711 

0.002715 

0.003031 

0.010604 

0.00 

1  micron 

0.002113 

0.000901 

0.002536 

0.004896 

1.73 

0.5  microns 

0.001604 

0.001573 

0.000706 

0.001017 

1.93 

Omicrons 

0 

0 

0 

0 

2 

Distance  Along  Asperity  (um) 


200mps 

*500mps 

lOOOmps 

•1500mps 


Figure  44:  Single  Asperity  formation  using  SMS  Failure  Criteria 


Table  13:  2um  Wear  per  Unit  Width  Data  for  Single  Asperity  Calculations  using  the  Johnson  and 

Cook  Failure  Criterion 


2um  Johnson  and  Cook  Failure 
Criteria 

Wear  rate 

oer  unit  widt 

h  (W^ 

200 

500 

1000 

1500 

Distance 

2  microns 

0.000907 

0.000771 

0.000940 

0.002536 

0.00 

1  micron 

0.000049 

0.000095 

0.000000 

0.000498 

1.73 

0.5  microns 

0.000000 

0.000043 

0.000036 

0.000069 

1.93 

0  microns 

0 

0 

0 

0 

2 
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Figure  45:  Single  Asperity  formation  using  JC  Failure  Criteria 

The  single  asperity  rates  for  each  of  these  velocities  was  obtained  by  finding  the 

average  area  under  the  curve  and  multiplying  by  2,  as  established  in  equation  4.3.  To  find 
the  correction  factor,  one  divides  the  single  asperity  rate  found  by  the  wear  per  unit  width 
of  the  radius  of  the  asperity.  Finally  by  taking  the  average  of  all  of  the  velocities  single 
asperity  rates,  one  arrives  at  a  value  that  can  be  used  generally  for  any  velocity.  Those 
results  are  seen  below  in  Tables  14  and  15. 

Table  14:  Spherical  Correction  Factor  Calculations  for  a  2um  asperity  using  the  Strain  at  Max  Stress 

Failure  Criterion 


Velocity 

Single  asperity  wear  rate 

(Wsa) 

Wuw  (2um) 

Ratio 

(Wsa/Wuw) 

m/s 

mm3/mm 

mm2/mm 

mm 

200 

9.203E-03 

2.71  E-03 

3.39E+00 

500 

6.862E-03 

2.72E-03 

2.53E+00 

1000 

1 .033E-02 

3.03E-03 

3.41  E+00 

1500 

2.809E-02 

1.06E-02 

2.65E+00 

Average: 

3.00E+00 
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Table  15:  Spherical  Correction  Factor  Calculations  for  a  2um  asperity  using  the  Johnson  and  Cook 

Failure  Criterion 


Velocity 

Single  asperity  wear  rate 

(Wsa) 

Wuw  (2um) 

Ratio 

(Wsa/Wuw) 

m/s 

mm3/mm 

mm2/mm 

mm 

200 

1.666E-03 

9.07E-04 

1.84E+00 

500 

1.530E-03 

7.71  E-04 

1.98E+00 

1000 

1.637E-03 

9.40E-04 

1.74E+00 

1500 

5.372E-03 

2.54E-03 

2.12E+00 

Average: 

1.92E+00 

4.5.2  Spherical  Correction  Factor  Comparison 

This  was  done  for  both  a  6um  asperity  and  a  2um  asperity.  The  6um  asperity  was 
originally  done  by  Chad  Hale  and  is  useful  for  comparing  to  previous  work.  Hale 
assumed  that  the  correction  factor  was  linearly  scalable  based  on  the  asperity  size.  The 
2um  correction  factor  was  of  interest  to  compare  to  Buentello’s  ABAQUS  simulation, 
where  a  2um  asperity  was  being  used  for  calculations.  All  of  the  calculations  performed 
to  arrive  at  the  final  numbers  are  included  in  Appendix  D.  Here  are  the  comparative 
results: 


Table  16:  Comparative  Spherical  Correction  Factors 


Source 

ABAQUS  [3] 

CTH 

CTH 

Failure  Criteria 

Johnson  and  Cook 

Strain  at  Max  Stress 

Johnson  and  Cook 

6um  Asperity 

Radius 

8.29 

10.1 

9.58 

2um  Asperity 

Radius 

2.76 

3.00 

1.92 

The  new  CTH  values  are  higher  than  the  ABAQUS  factor  the  6  micron  correction 
factor,  and  are  above  and  below  the  ABAQUS  values  for  the  2  micron  asperity.  The  CTH 
numbers  do  seem  reasonable  when  compared  to  the  ABAQUS  results.  This  is  particularly 
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encouraging  as  it  confirms  that  the  two  programs,  despite  their  differences,  can  produce 
competitive  results.  However  CTH’s  results  do  not  confirm  that  the  correction  factor  can 
be  scaled  based  on  asperity  size.  Now  that  the  spherical  correction  factor  has  been 
calculated,  wear  volumes  can  be  found. 


4.6  Wear  Results  and  Comparison 

For  the  results  section,  the  failure  criteria  will  be  abbreviated  to  SMS  (strain  at 
max  stress),  VM  (von  Mises),  and  JCF  (Johnson  and  Cook  Fracture)  as  a  space  saving 
measure.  The  wear  volume  rate  in  mm3/mm  is  found  using  the  following  equation: 


WVR 


WuwSCF 

Ds 


4.6 


Previous  work  by  Hale  and  Huber  used  6um  asperities,  so  for  the  moment  6um 
wear  results  will  only  be  used.  However,  now  there  are  different  failure  criteria  to  throw 
into  the  mix,  so  those  will  also  be  compared.  Below  are  the  results  found  in  this  study: 


Table  17:  Wear  Volume  Rates  using  Various  Failure  Criteria  and  SCF 


Wear  Volume  Rate  (mmA3/mm) 

ABAQUS  correction  Factor 
(8.29*10A-3) 

SMS  correction  Factor 
(10.1*10A-3) 

JCF  Correction  Factor 
(9.58*10A-3) 

SMS 

VM 

JCF 

SMS 

VM 

JCF 

SMS 

VM 

JCF 

200 

9.32E- 

05 

5.14E- 

05 

2.45E- 

05 

1.13E- 

04 

6.24E- 

05 

2.97E- 

05 

1.08E- 

04 

5.95E- 

05 

2.83E- 

05 

500 

6.16E- 

05 

5.49E- 

05 

3.08E- 

05 

7.48E- 

05 

6.66E- 

05 

3.74E- 

05 

7.13E- 

05 

6.35E- 

05 

3.57E- 

05 

1000 

8.72E- 

05 

2.78E- 

04 

4.02E- 

05 

1.06E- 

04 

3.37E- 

04 

4.87E- 

05 

1.01E- 

04 

3.22E- 

04 

4.64E- 

05 

1500 

9.99E- 

05 

3.29E- 

04 

4.65E- 

05 

1.21E- 

04 

3.99E- 

04 

5.64E- 

05 

1.16E- 

04 

3.80E- 

04 

5.37E- 

05 

Graphically,  they  look  something  like  this: 
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SCF=ABAQUS,  SMS 
SCF=ABAQUS,  VM 
SCF=ABAQUS/  JCF 
SCF=SMS,  SMS 
SCF=SMS,  VMS 
SCF=SMS,  JCF 
SCF=JCF,  SMS 


- SCF=JCF,  VM 

SMC=JCF,  JCF 


Figure  46:  Graphical  Wear  Volume  Comparison 

To  unpack  this  graph  a  bit,  there  seems  to  be  two  major  groups,  an  upper  wear 

bound  and  a  lower  wear  bound.  The  three  grouped  lines  at  the  very  bottom  of  the  graph 
are  wear  rates  using  the  Johnson  and  Cook  failure  criterion.  The  group  of  three  lines 
above  the  Johnson  and  Cook  failure  criterion  values  are  the  wear  rates  using  the  strain  at 
max  stress  failure  criteria.  Finally,  the  three  lines  grouped  at  the  top  of  the  graph  are  all 
the  wear  rates  generated  using  the  von  Mises  failure  criterion.  Applying  the  spherical 
correction  factors  provides  an  upper  bound  and  a  lower  bound  of  possible  wear  rates, 
with  the  strain  at  max  stress  correction  factor  providing  the  upper  bound  value  and  the 
Johnson  and  Cook  correction  factor  providing  the  lower  bound  value. 

The  Johnson  and  Cook  wear  rates  are  all  lower  than  the  von  Mises  and  strain  at 
max  stress  wear  rates.  This  is  because  of  the  feature  of  the  Johnson  and  Cook  failure 
criterion  to  allow  the  strength  of  damaged  materials  to  be  set  to  zero  during  the 
simulation.  Removing  the  material  strength  lessens  the  pressure  in  the  slipper  and  reduces 
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the  overall  wear.  When  material  is  not  allowed  to  leave,  pressure  accumulates  more  and 
the  wear  rate  is  increased. 

How  do  these  results  stack  up  against  previous  work  results?  Here  is  a  visual 
comparison  of  this  research’s  work  (prefaced  with  LW  in  the  legend,  the  author’s  initials) 
to  Huber,  Hale,  and  Wolfson’s  data: 


E 

E 


ro 

cc 


ro 

cu 


<u 

E 

O 

> 


Dave  Huber  SMS 
Dave  Huber  VM 
Chad  HaleJCF 
Wolfson  Data 
LW  SMS 
LW  VM 
LWJCF 


Figure  47:  Comparison  of  Wear  Volume  Rates 
For  this  comparison,  only  wear  rates  using  the  ABAQUS  spherical  correction 

factor  were  used  in  an  attempt  to  reduce  clutter  on  the  graph  and  show  the  median  of 

possible  wear  volume  rate  results.  The  new  strain  at  max  stress  and  von  Mises  wear 

volume  rates  line  up  very  closely  with  the  work  done  by  Huber.  The  new  values  differ 

from  previous  work  with  the  CTH  model  because  the  simulation  in  this  research  was  run 

for  a  longer  period  of  time.  In  the  slower  velocity  regimes,  Huber’s,  Hale’s,  and  the  new 
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strain  at  max  stress  and  von  Mises  values  seem  to  be  rather  close  to  the  Wolfson 


experimental  data.  The  Johnson  and  Cook  failure  criterion  seems  to  under-predict  the 
wear  here.  This  may  be  an  indication  that  friction  wear  is  more  relevant  in  the  200mps 
velocity  regime.  However,  as  the  velocity  increases,  the  von  Mises  wear  increases  much 
more  compared  to  the  other  strain  rates.  There  is  also  a  significant  difference  in  the  strain 
at  max  stress  failure  criterion  results  from  previous  research.  This  clearly  indicated  that 
finding  when  the  maximum  damage  area  occurs  will  change  the  results.  The  difference  is 
not  as  dramatic  in  between  the  von  Mises  criteria  because  the  von  Mises  damage  area 
reaches  an  asymptote  earlier  than  any  of  the  other  criterion  being  used. 

Finally  it  is  of  interest  to  try  to  compare  with  Rodolfo  Buentello’s  3D  ABAQUS 
model.  This  study’s  results  can  compare  to  his  calculations  at  500mps.  In  his  work  he 
found: 


Table  18:  Buentello  3D  ABAQUS  Damage  Volume  Calculations 


Stage 

Ave. 

Vel 

(m/s) 

WuA 

(mmVmm) 

Width 

(mm) 

Ave. 

Worn 

length 

Sliding 

Distance 

% 

Contact 

Predicted 

Volume 

Damage 

(mm3) 

Measured 

Volume 

Damage 

(min' 

2 

500 

1.81E-6 

101 

35 

1016000 

0.30 

1993.9 

1684.2 

Where  the  volume  of  damage  was  found  using  the  following  equation: 

n 

vd  =  ^  WUAnACnDsn%Contn  3.2 

n= 1 

To  be  able  to  compare  to  these  numbers,  the  wear  volume  per  unit  area  is  needed. 
Using  an  equation  developed  by  Buentello,  the  following  equation  is  used: 
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4.7 


WuA  = 


Ad  SCF  x  10 
Measured  Area 


Where  Ad  is  the  damage  area,  the  measured  area  is  the  0.16mm~  and  10  is  the 


number  of  asperities  found  in  the  0.16  area 


The  results  using  each  of  the  spherical  correction  factors  are  shown  below: 


Table  19:  Wear  Area  Rate  Calculations 


Fracture 

Criteria 

mm3/mm2 

mm3 

3D ABAQUS 

JCF 

1.81E-06 

2.90E-08 

SCF  =  SMS 

SMS 

5.15E-06 

8.25E-08 

VM 

2.70E-06 

4.33E-08 

JCF 

1.46E-06 

2.34E-08 

SCF  =  JCF 

SMS 

3.30E-06 

5.29E-08 

VM 

1.73E-06 

2.77E-08 

JCF 

9.39E-07 

1.50E-08 

SCF  =  ABAQUS 

SMS 

4.76E-06 

7.61E-08 

VM 

2.49E-06 

3.99E-08 

JCF 

1.35E-06 

2.16E-08 

Interestingly  enough,  the  c 


iosest  CTH  results  are  the  ones  using  the  von  Mises 


failure  criterion.  They  actually  are  extremely  close. 


It  is  also  possible  to  attempt  to  match  the  volume  of  wear  that  Buentello 


predicted.  Using  Equation  3.2  and  using  the  same  sliding  distance,  %contact,  and  area  of 


contact,  the  following  results  are  generated: 


Table  20:  Comparative  Predicted  Volume  Damage 


Spherical 

Correction 

Factor  used 

Fracture 

Criteria  used 

Wear  per  unit 

area 

(mm3/mm2) 

Predicted 

Volume 

Damage 

Measured 

Volume 

Damage 

% 

Difference 

ABAQUS 

SMS 

4.76E-06 

5127.5 

1684.2 

204.5 

VM 

2.49E-06 

2687.2 

1684.2 

59.6 

JCF 

1.35E-06 

1456.6 

1684.2 

13.5 

SMS 

SMS 

5.1581E-06 

5557.7 

1684.2 

230.0 

VM 

2.7032E-06 

2912.6 

1684.2 

72.9 
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JCF 

1.46526E-06 

1578.8 

1684.2 

6.3 

JCF 

SMS 

3.30623E-06 

3562.4 

1684.2 

111.5 

VM 

1.73269E-06 

1866.9 

1684.2 

10.84 

JCF 

9.39201E-07 

1011.9 

1684.2 

39.9 

3D  simulation 

JCF 

1.81E-06 

1950.2 

1684.2 

15.8 

Admittedly,  most  of  the  results  are  not  exactly  spot-on.  However,  for  the  two  of 
the  Johnson  and  Cook  wear  rates,  the  results  are  actually  better  than  those  predicted  with 
the  3D  ABAQUS  simulation. 


The  comparison  seems  to  indicate  that  some  failure  criteria  may  be  more 
appropriate  for  modeling  certain  velocity  ranges  than  others.  Certainly  in  the  comparative 
data,  the  Johnson  and  Cook  data  was  predicting  values  much  lower  than  the  experimental 
data  at  lower  values,  but  the  SMS  and  von  Mises  data  was  much  closer.  It  seems  that 
keeping  the  material  from  being  removed  in  the  simulation  can  actually  compensate  for 
the  lack  of  frictional  representation  at  lower  velocities.  However,  using  a  combination  of 
failure  criteria  is  rather  cumbersome,  and  not  exactly  the  streamlined  simple  solution  that 
is  sought.  Also  this  assumption  is  based  on  the  3D  Abaqus  data  behaves  like  the  Wolfson 
experimental  data,  which  is  a  relationship  that  has  not  officially  been  established.  Thus 
this  postulation  cannot  be  confirmed  without  additional  data  to  back  it  up. 


4.7  Wear  Results  Summary 

To  summarize  what  was  done  in  this  section,  the  background  for  the  wear  per  unit 
width,  single  asperity  wear  rate,  spherical  correction  factor,  and  wear  volume  were  all 
reviewed,  and  results  of  these  calculations  were  presented.  It  was  show  that  during  the 
simulation,  a  maximum  damage  area  will  occur,  but  not  all  of  the  maximum  damage 
areas  have  been  found  yet  using  CTH.  Following  this,  single  asperity  rates  were 
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calculated,  which  showed  that  wear  of  a  single  asperity  increases  with  increasing  velocity 
and  asperity  size.  The  spherical  correction  factor  was  calculated,  and  the  results  found  by 
Hale’s  plane  strain  model  in  ABAQUS  within  the  range  of  results  produced  by  the  CTH 
Johnson  and  Cook  failure  criterion.  The  strain  at  max  stress  correction  factor  produced 
higher  values,  as  was  expected  since  there  is  no  material  removal  with  that  criterion. 

Wear  volume  rates  were  calculated  using  the  new  spherical  correction  factors  and 
compared  to  past  work  done  by  Huber,  Hale,  and  Wolfson.  The  wear  volume  rates 
calculated  using  the  old  criteria  matched  up  with  previous  work  fairly  well,  and  the 
Johnson  and  Cook  criterion  produced  lower  wear  volume  rates,  which  were  closer  in 
value  to  the  Wolfson  experimental  data  for  higher  velocities.  At  lower  velocities,  the 
strain  at  max  stress  and  von  Mises  wear  volume  rates  appeared  to  be  closer  to  the 
experimental  Wolfson  Data  than  the  Johnson  and  Cook  wear  volume  rates.  Finally,  wear 
area  rates  and  wear  volume  was  calculated  and  compared  to  a  3D  ABAQUS  simulation 
results.  The  Johnson  and  Cook  failure  criterion  seemed  to  produce  the  best  results,  but 
produced  smaller  results  when  used  in  combination  with  the  Johnson  and  Cook  spherical 
correction  factor. 
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V.  Conclusions 


5.1  Research  Summary 

This  project  set  about  attempting  to  demonstrate  the  validity  of  a  2D  plane  strain 
hydrocode  model  to  model  the  wear  of  a  slipper  in  a  high  velocity  environment.  The  CTH 
model  was  advanced  by  engaging  the  Johnson  and  Cook  failure  criterion,  by  re¬ 
calculating  spherical  correction  factors  previously  developed  by  Chad  Hale,  and  by  using 
a  wear  data  collected  from  longer  simulation  runs.  The  study  shows  that  the  wear  results 
are  comparable  to  those  calculated  by  ABAQUS. 

5.2  Conclusions 

This  study  has  allowed  us  to  reach  several  conclusions.  It  was  our  goal  to 
determine  if  using  a  2D  plane  strain  scenario  with  the  Johnson  and  Cook  fracture 
criterion  would  improve  previous  models  predictions  of  wear.  Happily  it  can  be  said  that 
that  the  results  showed  an  overall  improvement  of  wear  results  closer  to  those  seen  in  an 
ABAQUS  3D  model  over  other  failure  criteria  used  previously. 

One  of  the  questions  that  this  model  raises  is  if  the  asperity  collision  is  the  driving 
factor  of  the  wear,  since  in  this  CTH  model  friction  cannot  contribute  anything.  Based  on 
the  results,  asperity  collisions  can  account  for  the  majority  of  the  wear,  as  it  predicts  wear 
closer  than  those  seen  in  experimental  data.  More  data  from  ABAQUS  is  needed  to 
confirm  this  when  comparing  simulation  models. 

There  had  been  concern  that  a  2D  plane  strain  model  might  be  over-simplified 
and  not  be  able  to  capture  enough  of  the  wear  event  to  make  it  a  competitive  model. 
However  since  the  results  are  in  the  range  as  those  produced  in  3D  simulations  and 
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compared  to  data  from  HHSTT,  it  seems  the  simplified  model  will  work  well  given  the 
right  failure  criterion  is  in  place. 

One  theory  that  was  confirmed  was  that  engaging  the  Johnson  and  Cook  failure 
criterion  in  CTH  would  lead  to  lower  predicted  wear  than  previously  found  using  other 
failure  criteria.  This  in  fact  turned  out  to  be  true  and  was  reflected  in  the  results  in  that  the 
Johnson  and  Cook  wear  was  lower  than  the  Wolfson  data  and  ABAQUS  data. 

Some  would  worry  about  how  far  one  can  extend  a  single  asperity  collision  to  an 
overall  wear  event,  but  our  results  coincide  with  the  previous  research  that  concludes  that 
this  is  reasonable  approach  that  produces  realistic  results. 

One  thing  that  can  be  noted  is  that  CTH,  while  fantastic  at  speeds  greater  than 
500mps,  seems  to  misbehave  at  lower  velocities.  It  is  very  promising  for  research  that 
will  want  to  go  to  speeds  greater  than  those  in  this  research,  but  for  lower  velocities,  it 
may  be  best  that  another  program  that  can  model  friction  better  be  used. 

When  determining  the  spherical  correction  factor,  it  was  found  that  CTH  and 
ABAQUS  both  created  similar  results.  This  is  fantastic  because  it  shows  how  compatible 
the  two  modeling  programs  are,  and  further  confirms  that  they  both  should  predict  similar 
values  despite  using  different  solving  methods.  Using  the  CTH-determined  spherical 
correction  factors  provided  an  upper  and  lower  bound  for  wear  rate  calculations. 
Additionally,  it  was  found  that  the  idea  that  the  spherical  correction  factor  can  be  scaled 
based  on  asperity  size  may  not  actually  be  true. 

In  most  cases,  the  wear  rate  increased  with  increasing  velocities. 
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5.1  Future  Work 


As  with  most  research,  as  results  come  in,  more  questions  arise.  In  the  future  it 
would  be  beneficial  to  further  investigate  the  differences  in  the  ABAQUS  3D  model  and 
the  CTH  model  for  additional  velocities  and  asperity  sizes. 

For  the  failure  criterion  damage  areas,  simulations  will  need  to  be  run  longer  to 
find  when  the  maximum  damage  occurs.  Once  that  is  done,  it  would  be  interesting  to  take 
the  data  and  see  if  the  maximum  damage  occurred  at  predictable  distances.  This  would  be 
helpful  so  that  one  would  not  have  to  guess  at  how  long  a  simulation  needs  to  run  to 
obtain  the  maximum  damage. 

While  the  spherical  correction  factor  was  validated  between  plane  strain  models 
in  CTH  and  ABAQUS,  it  would  be  interesting  to  do  more  work  in  the  3D  coordinate 
system.  CTH  has  3D  capabilities,  and  a  3D  model  is  already  in  development.  It  would  be 
very  interesting  to  finish  developing  this  model  and  make  comparisons  between  plane 
strain  simulation  and  the  3D  ABAQUS  simulation. 

Additionally,  there  are  several  modifications  that  might  be  made  to  the  model 
itself.  There  is  work  going  on  investigating  the  true  initial  temperature  of  the  materials  in 
the  model,  so  it  would  be  of  interest  implement  a  more  accurate  number  there.  Another 
possible  modification  would  be  to  move  away  from  the  Johnson  and  Cook  Constitutive 
Model  and  investigate  other  potential  constitutive  models  and  make  comparisons.  An 
additional  modification  would  be  implementing  the  Mie-Griineisen  Equation  of  state  in 
the  CTH  model  in  place  of  the  SESAME  tables  to  see  if  results  closer  match  the 
ABAQUS  predictions. 
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Appendix  A  -  CTH  Input  File 


Below  is  the  input  file  that  was  imported  into  CTH  for  this  research.  Parameters  that  were 
changed  were  the  velocity,  simulation  times,  and  size  of  asperity,  along  with  turning  off 
the  Johnson  and  Cook  fracture  criterion  by  commenting  out  its  parameters.  It  has  been 
updated  to  be  compatible  with  CTH  version  10.2. 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'kit'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'kit'k'k'k'k 

*eor*  cthin 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

*  cthin  input  with  Spymaster  graphics  for  slipper  wear  simulation 

*  filename:  slipperwear . in 

*  1.  File  modified  by  Steve  Meador  (MS-10M) 

*  2.  File  converted  to  CTH  v8 . 1  by  Maj  Chad  Hale,  PhD-09S,  Aug 
2008 

*  3.  new  format  based  on  CTH  Course  (4-7  Aug  08)  in  Albuquerque, 
NM 

*  4.  modifies  Cameron's  393  m/s.  No  Coating,  Asperity,  T=297 
input  file 

9c 

9c 

*  |  - >  | 

*  I  I  I 

*  I  v  / 

9r  _ 

9c 

9c 

*  vx=varies,  vy=-l  m/s  V300  Steel  Slider,  1080  Steel  Rail,  No 
Atm. 

*  No  Slide  line.  mix=l  frac=l  Rounded  corner. 

*  Added  mass  on  top  to  simulate  sled  mass 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

*  title  record  set 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

Horizontal  Velocity  =  500  m/s,  Vertical  Velocity  =  -0.50  m/s 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

*  control  input  set 

9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c9c 

control 
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mmp3  *  enable  multiple  material  temperatures 

and  pressures  in  each  cell 

tstop  =  4.40e-9  *  stopping  criteria  for  time  level  -  this 

is  total  simulation  time 

nscycle  =  100000  *  maximum  number  of  cycles  to  be  run 

*  rdumpf  =  3600.  *  time  for  back-ups  of  restart  file 

updates 

tbad  =  le30  *  maximum  number  of  thermodynamics 

warnings 

*  dtcourant  =0.6  *  Courant  condition  multiplier 

ygravity  =  -980  *  Acceleration  due  to  gravity  =  -9.80 

m/sA2 

endcontrol 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'kic 

*  mesh  input  set 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

*  geom=2DR (rectangular  x,y) 

*  geom=2DC ( cylindrical  x=radius,  y=axis) 

*  geom=3DR (rectangular  x,y,z) 

*  type=e  (Eulerian)  -  now  the  default  (CTHv8.1) 

*  x#=coordinate  range  for  plot 

*  y#=coordinate  range  for  plot 

*  dxf=width  of  first  cell  in  the  region 

*  dxl=width  of  last  cell  in  the  region 

*  n=number  of  cells  added  in  this  region 

*  w=total  width  of  this  region  in  centimeters 

*  r=ratio  of  adjacent  cell  widths 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

mesh 

block  1  geom=2dr 
Eulerian  mesh 

xO  =  0.0000 

xl  w  =  850e-4  dxf 
endx 

yO  =  0.0000 

yl  w  =  850e-4  dyf  =  1 . Oe-4  dyl  =  1 . Oe-4 
endy 

endblock 
endmesh 

'k'k'k'k'k'k'kit'k'k'k'k'k'k'k'kit'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'kit'k'k'k'k 

*  EOS  input  set 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

eos 

materiall  ses  iron  *  1080  steel  rail 

material2  ses  steel_v300  *  VascoMax  300  slipper 

endeos 


*  coordinates  for  2D  rectangular 

=  1 . Oe-4  dxl  =  1 . Oe-4 
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ic'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 


*  elastic-plastic  input  set 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

epdata 

vpsave  *  cell  yield  stress  and  plastic  strain  rate  data 

is  saved 

lstrain  *  compute  and  save  Lagrangian  strain  tensor 

components 

mix  =3  *  volume  averaged  yield  strength  normalized  by 

sum  of  volume  fractions 


matep  =1  *  1080 

JO  USER 
A JO  0 . 7el0 
BJO  3 . 6el0 
CJO  0.17 
MJO  0.25 
NJO  0.6 
T JO  0.14391 
poisson  0.27 


Steel  rail 

*  A 

*  B 

*  C 

*  m 

*  n 

*  Melting  temperature 


matep  =2  *  VascoMax  300 


JO  USER 

A JO  =  2  .  lelO 

~k 

A 

BJO  =  0 . 124el0 

~k 

B 

CJO  =0.03 

~k 

C 

MJO  =0.8 

~k 

m 

NJO  =  0 .3737  *  n 

T JO  =  0.145202 

~k 

Melting 

poisson  0.283 

JFRAC  USER 

JFD1  =  -0.09 

~k 

D1 

JFD2  =  0.27 

D2 

JFD3  =0.48 

~k 

D3 

JFD4  =  0.014 

D4 

~k 

JFD5  =  3.870 

~k 

D5 

JFTM  =  0 .145202 

Melting 

JFPF0  =  -7 . 45el0 

=PFRAC2 

slipper 


temperature 


temperature 


SLIDE  1  2 
endepdata 
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'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 


*  diatom  input  set 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

diatom 
block  1 

package  '1080  steel  rail' 
material  1 
numsub  100 

temperature  =  2.55935e-2  *  eV  =  74.93F  =  297  K 
velocity  0.0,  0.0 
insert  box 
pi  0  0 

p2  850e-4  200e-4 
endinsert 
delete  circle 
center  695e-4  200e-4 
radius  2e-4 
enddelete 
insert  circle 
center  695e-4  200e-4 
radius  2e-4 
endinsert 
endpackage 

package  'slipper' 
material  2 
numsub  100 

temperature  =  0.0184558 
velocity  =  500e2,  -0.50e2 
insert  box 

pi  0.0  200e-4 
p2  694e-4  325e-4 
endinsert 
delete  box 
pi  692e-4  200e-4 
p2  694e-4  202e-4 
enddelete 
delete  circle 
center  692e-4  202e-4 
radius  2e-4 
enddelete 
insert  circle 
center  692e-4  202e-4 
radius  2e-4 
endinsert 
endpackage 

endblock 

enddiatom 
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'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 


*  tracer  input  set 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

tracer 

add  0.06755,  0.01905  to  0.07115,  0.01905  n=37 

add  0.06755,  0.01915  to  0.07115,  0.01915  n=37 

add  0.06755,  0.01925  to  0.07115,  0.01925  n=37 

add  0.06755,  0.01935  to  0.07115,  0.01935  n=37 

add  0.06755,  0.01945  to  0.07115,  0.01945  n=37 

add  0.06755,  0.01955  to  0.07115,  0.01955  n=37 

add  0.06755,  0.01965  to  0.07115,  0.01965  n=37 

add  0.06755,  0.01975  to  0.07115,  0.01975  n=37 

add  0.06755,  0.01985  to  0.07115,  0.01985  n=37 

add  0.06755,  0.01995  to  0.07115,  0.01995  n=37 

add  0.06755,  0.02005  to  0.07115,  0.02005  n=37 

add  0.06755,  0.02015  to  0.07115,  0.02015  n=37 

add  0.06755,  0.02025  to  0.07115,  0.02025  n=37 

add  0.06755,  0.02035  to  0.07115,  0.02035  n=37 

add  0.06755,  0.02045  to  0.07115,  0.02045  n=37 

add  0.06755,  0.02055  to  0.07115,  0.02055  n=37 

add  0.06755,  0.02065  to  0.07115,  0.02065  n=37 

add  0.06755,  0.02075  to  0.07115,  0.02075  n=37 

add  0.06755,  0.02085  to  0.07115,  0.02085  n=37 

add  0.06755,  0.02095  to  0.07115,  0.02095  n=37 

add  0.06755,  0.02105  to  0.07115,  0.02105  n=37 

add  0.06755,  0.02115  to  0.07115,  0.02115  n=37 

add  0.06755,  0.02125  to  0.07115,  0.02125  n=37 

add  0.06755,  0.02135  to  0.07115,  0.02135  n=37 

add  0.06755,  0.02145  to  0.07115,  0.02145  n=37 

add  0.06755,  0.02155  to  0.07115,  0.02155  n=37 

add  0.06755,  0.02165  to  0.07115,  0.02165  n=37 

add  0.06755,  0.02175  to  0.07115,  0.02175  n=37 

add  0.06755,  0.02185  to  0.07115,  0.02185  n=37 

add  0.06755,  0.02195  to  0.07115,  0.02195  n=37 

add  0.06755,  0.02205  to  0.07115,  0.02205  n=37 

add  0.06755,  0.02215  to  0.07115,  0.02215  n=37 

add  0.06755,  0.02225  to  0.07115,  0.02225  n=37 

add  0.06755,  0.02235  to  0.07115,  0.02235  n=37 

add  0.06755,  0.02245  to  0.07115,  0.02245  n=37 

endtracer 


icicicicicicic'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

*  convection  control  input  set 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

Convct  *  enable  convection  of  internal  energy 

convection  =1  *  use  slope  of  internal  energy  and  mass 

density,  discard  KE  residual 

interface  =  smyra  *  scheme  for  interface  tracker 
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endconvct 


icit'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

*  fracture  input  set 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

Fracts  *  enable  fracture  data  (dynes/cmA2) 

pressure 

pfracl  =  -2.0el0 
pfrac2  =  -7.45el0 
endf racts 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

*  edits  input  set 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 


edit 

exact 

shortta 

time  =  0.0 
ends 
longt 

time  =  0 . OeO 
endl 
plott 

time  0 . Oe-6 
endp 
histt 

time  0 . Oe-6 
htracer  all 
endhistt 
ende 


*  short  edits  based  on  time 
dt  =  4 . 40e-ll 

*  long  edits  based  on  time 
,  dt  =  4 . 40e-ll 

*  plot  dumps  based  on  time 
dtfrequency  4.40e-ll 

*  tracer  history  based  on  time 
dtfrequency  4.40e-ll 


icicic'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

*  boundary  condition  input  set 

'k'k'k'k'k'k'k'kit'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

*  0=symmetry 

*  l=sound  speed  based  absorbing 

*  2=extrapolated  pressure  with  no  mass  allowed  to  enter 

*  3=extrapolated  pressure  but  mass  is  allowed  to  enter 

'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 


boundary 
bhydro 
conditions 
block  1 

bxbot  =  1  ,  bxtop  =  2 
bybot  =  1  ,  bytop  =  2 
endb 
endh 
endb 


*  enable  boundary  condition  data 

*  enable  hydrodynamic  boundary 
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maxdt  *  maximum  allowable  time  step  in  mesh 

time  =  0.0  dt  =  4.40e-ll 
endm 


*  CSH:  Attempt  to  get  data  for  Spymaster 


spy 

PlotTime ( 0 . 0,  4.40e-ll); 
SaveTime ( 0 . 0,  4.40e-ll); 


Save ( "VOID, VOLM, M, P, XXDEV, YYDEV, XYDEV,  VX,  VY,  T,  TK,  PM,  TM,  YLD,  Q2  ,  J2P 


define  main ( ) 


%  pprintf ( "  PLOT:  Cycle=%d,  Time=%e\n ",  CYCLE, TIME ) ; 

%  XLimits ( 400e-4, 725e-4) ; 

%  YLimits ( 175e-4, 300e-4)  ; 

%  Image ( "Materials ")  ; 

%  Window ( 0 ,  0 ,  0 . 75 ,  1 )  ; 

%  Label ( sprintf ( "Materials  at  %6.2e  seconds",  TIME)); 

%  Plot2DMats (0.3) ; 

%  ULabel ( "Test :  (cm)"); 

%  Draw2DMesh ( ) ;  %  toggle  on/off  mesh 

%  MatColors (RED, GREEN, YELLOW, NO_COLOR) ; 

%  MatNames ( "Epoxy  Coating", "1080  Steel  Rail " , "VascoMax  300 
Slipper", " ") ; 

%  Dr awMat Legend ("",0.71,0.2,0.99,0.9); 

%  Endlmage; 


XLimits (650e-4,  750e-4)  ; 

YLimits (150e-4,250e-4)  ; 

Image ( "VonMisesStress" )  ; 

Window (0,0,0.75,1); 

ColorMapRange (0,4000); 

ColorMapClipping (OFF, OFF) ; 

Label ( sprintf ( "von  Mises  Stress  at  %6 . 2e  seconds",  TIME)); 
Plot2D ( " J2P" ) ; 

Dr aw2DMat Contour ; 

DrawColorMap ( " vonMises  Stress  (MPa)",  0 . 7,  0 . 4,  0 . 9,  0 . 9) ; 
Endlmage ; 

XLimits (650e-4, 750e-4) ; 

YLimits (150e-4,250e-4)  ; 

Image ( "PlasticStrainRate"  )  ; 

Window (0,0,0.75,1); 
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%  ColorMapRange ( le6 , lel5, LOG_MAP) ; 

%  ColorMapClipping (OFF,  OFF)  ; 

%  Label ( sprintf ( "Plastic  Strain  Rate  at  %6.2e  seconds", 
TIME)  )  ; 

%  Plot2D ( "PSR"  )  ; 

%  Draw2DMatContour ; 

%  DrawColorMap ( "Plastic  Strain  Rate  (1/sec)", 

0 . 7, 0. 4, 0. 9, 0. 9)  ; 

%  Endlmage; 

XLimits (685e-4,  715e-4)  ; 

YLimits (190e-4,  215e-4)  ; 

Image ( "Materials_small " ) ; 

Window (0,0,0.75,1); 

Label ( sprintf ( "Materials  at  %6.2e  seconds",  TIME)); 
Plot2DMats (0.3); 

Label (  "Test  Label:  Distance  (cm)"  ); 

%  Draw2DMesh ( ) ;  %  toggle  on/off  mesh 

MatColors (GREEN, YELLOW) ; 

MatNames ( " 10 80  Steel  Rail " , "VascoMax  300  Slipper"); 

Dr awMatLegend ("",0.71,0.2,0.99,0.9)  ; 

Endlmage; 

XLimits (650e-4, 750e-4) ; 

YLimits (150e-4,  250e-4)  ; 

Image ( "Pressure"  )  ; 

Window (0,0,0.75,1); 

ColorMapRange ( le6 , 2ell , LOG_MAP ) ; 

ColorMapClipping (OFF, OFF) ; 

Label ( sprintf ( "Pressure  at  %6.2e  seconds",  TIME)); 
Plot2D ( "P" ) ; 

Draw2DMatContour ; 

DrawColorMap ( "Pressure  (dyne/cmA2 ) " ,  0 . 7, 0 . 4, 0 . 9, 0 . 9) ; 
Endlmage; 


SaveHis ( "POSITION, Q2 , PSR, VOLM+2 ,  J2P,  YLD" )  ; 

%SaveHis ( "POSITION, Q2 , PSR, VOLM+2 , J2P, DMG2 " ) ; 

% SaveHis ("POSITION, YLD, Q3 , PSR, VOLM+3 , P, XXDEV, YYDEV, J2P, DMG3 " ) ; 
%  SaveHis ( "POSITION, YLD, VOLM+3 , P , XXDEV, YYDEV, XYDEV" ) ; 

%  SaveHis ( "POSITION, Q3 , PSR, VOLM+3 " ) ; 

%  SaveHis ( "POSITION, VOLM+3, DMG3") ; 

SaveTracer (ALL) ; 

His Time (0,4.40e-ll)  ; 

define  spyhis_main ( ) 

{ 

HisLoad(l, "hscth") ; 

Label ("EFP  Velocity  (Tracer  1)"); 
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TPlot ( "VY . 1 " , 1 , AUTOSCALE ) ; 

} 

endspy 


Appendix  B.l  Strain  at  Max  Stress  MATLAB  Post-Processing  File 


Below  is  the  post-processing  file  used  to  calculate  the  wear  per  unit  with  per  time  for  the 
strain  at  max  stress  failure  criterion. 

%  Lauren  Wuertemberger 
%  AFIT 

%  7  Aug  2012 
%  Importing  Wear  data 

clear  all 
close  all 
clc 

velocity  =  500; 
meshSize  =  1 . 0e-4*l . Oe-4 ; 

%  Load  wear  data  file 

dataFile  =  '2um500mpsData-2.txt'; 

data  =  load (dataFile ) ; 

disp( 'Data  Imported. . . ' ) 

%  Start  grabbin '  stuff 

time  =  data ( : , 1 ) ; 

CycleNum  =  length (time) ; 

TracerPoints  =  ( size (data, 2 ) -1 ) / 7 ; 
xLoc  =  2; 
yLoc  =  3; 
yieldLoc  =  4; 

J2PLoc  =  5; 
vfLoc  =  6; 
srLoc  =  7;  %PSR 
sLoc  =8;  %Q 

for  iter  =  1 : TracerPoints 

xPoints ( : ,  iter )  =  data ( : , xLoc ) ; 
yPoints ( : , iter )  =  data (:, yLoc) ; 
yieldData ( : ,  iter )  =  data (:, yieldLoc ) ; 

J2PData ( : , iter )  =  data ( : , J2PLoc ) ; 
vf Data (:, iter )  =  data (:, vfLoc ) ; 
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srData ( : , iter )  =  data ( : , srLoc) ; 
strainData ( : , iter )  =  data ( : , sLoc ) ; 

xLoc  =  xLoc  +  7 ; 
yLoc  =  yLoc  +  7; 
yieldLoc  =  yieldLoc+7; 

J2PLoc  =  J2PLoc  +  7;  %% 
vfLoc  =  vfLoc  +  7; 
srLoc  =  srLoc  +7;  %% 
sLoc  =  sLoc  +  7; 

end 

%  EVALUATE  STRAIN  RATES  FOR  ZEROS 

for  r  =  1 : size ( srData, 1 ) 

for  c  =  1 : size ( srData, 2 ) 
if  srData (r, c) <. 002 

srData(r,c)  =  .002; 

end 

if  srData (r, c) >10el7 

disp ( ' Temp : ' ) , disp (temp) 
disp('H  Vel :'), disp (velocity ) 
disp('V  Vel :'), disp (vVel ) 
disp ( ' Row : ' ) , disp (r ) 
disp ( ' Col : ' ) , disp ( c ) 

disp ( 'Strain  Rate ' ) , disp ( srData (r , c ) ) 
error (' Strain  Rate  Out  of  Range') 

end 

end 

end 


%  EVALUATE  STRAIN  AT  MAX  STRESS  FAILURE  AREA 

failureSMS  =  zeros (CycleNum, TracerPoints ) ; 

f ailureSumSMS  =  zeros (CycleNum, 1 ) ; 

A  =  2 . 24700e-2 ; 

B  =  -5.5160e-2; 

C  =  6 . 04400e-3 ; 

f ailureCritSMS  =  A* ( srData . AB)  +  C; 

for  row=l : r 

for  col=l : c 

if  row>l  &&  failureSMS (row-1,  col) ==1 
failureSMS (row, col ) =1 ; 

end 

if  strainData (row,  col ) >=f ailureCritSMS (row, col ) 
failureSMS (row, col ) =1 ; 
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end 


end 

end 

failureSMS  =  failureSMS . *vfData;  %unitless 

for  iter  =  1 : length ( failureSumSMS ) 

f ailureSumSMS ( iter , 1 )  =  sum ( failureSMS ( iter ;  %unitless 

end 

damAreaSMS  =  f ailureSumSMS*meshSize ;  %mmA2 
for  i  =  2:CycleNum 

distanceSlid  =  velocity*time ( i ) *1000 ; 

WR_SMS ( i )  =  100*damAreaSMS (i) /distanceSlid; 

end 

SMS_Average  =  100*damAreaSMS (end) /distanceSlid 

%WR_SMS  =  100*damAreaSMS/distanceSlid;  %mmA2/mm 

disp(' Strain  at  Max  Stress  Failure  Mechanism  Evaluated...') 

%Calculate  Wear  at  each  interval 
plot (time,  WR_SMS) 
xlabel('Time  (s)') 

ylabel('Wear  Rate  per  Unit  Width  (mmA2/mm) ' ) 
title ( ' 2um  Asperity  Wear') 

figure 

plot (time, damAreaSMS ) 
xlabel ( ' Time  ( s ) ' ) 

ylabel('Area  of  Damaged  Elements (mmA2 )' ) 
title ( ' 2um  Asperity  Damage  Area') 


Appendix  B.2  Johnson-Cook  Failure  Criteria  Post-Processing  MATLAB  File 


Below  is  the  post-processing  file  used  to  calculate  the  wear  per  unit  with  per  time  for  the 
Johnson  and  Cook  failure  criterion. 


%  Lauren  Wuertemberger 
%  AFIT 

%  7  Aug  2012 
%  Importing  Wear  data 

clear  all 
close  all 


99 


clc 


velocity  =  1500; 
meshSize  =  1 . 0e-4*l . Oe-4; 

%  Load  wear  data  file 

dataFile  =  '2uml500mpsJFData-2.txt'; 

data  =  load (dataFile ) ; 

disp('Data  Imported...') 

%  Start  grabbin '  stuff 

time  =  data ( : ,  1 ) ; 

CycleNum  =  length (time) ; 

TracerPoints  =  ( size (data, 2 ) -1 ) /8 ; 
xLoc  =  2; 
yLoc  =  3; 

JFLoc  =  4; 
yieldLoc  =  5; 

J2PLoc  =  6; 
vfLoc  =  7; 
srLoc  =  8;  %PSR 
sLoc  =9;  %Q 

for  iter  =  1 : TracerPoints 

xPoints ( :  ,  iter )  =  data ( : , xLoc) ; 
yPoints ( : , iter )  =  data (:, yLoc ) ; 

JFData ( : , iter )  =  data (:, JFLoc) ; 
yieldData ( : , iter )  =  data (:, yieldLoc ) ; 
J2PData ( : , iter )  =  data ( : , J2PLoc ) ; 
vf Data (:,  iter )  =  data (:, vfLoc) ; 
srData ( : , iter )  =  data (:, srLoc ) ; 
strainData ( : , iter )  =  data ( : , sLoc ) ; 

xLoc  =  xLoc  +  8 ; 
yLoc  =  yLoc  +  8; 

JFLoc  =  JFLoc  +  8; 
yieldLoc  =  yieldLoc+8; 

J2PLoc  =  J2PLoc  +8;  %% 
vfLoc  =  vfLoc  +  8 ; 
srLoc  =  srLoc  +8;  %% 
sLoc  =  sLoc  +  8; 

end 

%%  EVALUATE  STRAIN  RATES  FOR  ZEROS 

O, 

O 

%  for  r  =  1 :  size ( srData, 1 ) 

%  for  c  =  1 :  size ( srData, 2 ) 

%  if  srData (r, c) <. 002 

%  srData(r,c)  =  .002; 

%  end 

o, 

o 

%  if  srData (r, c) >10el7 


100 


o\o  o\o  o\o 


o, 

o 

o, 

o 

o, 

o 

g, 

o 

g, 

o 

g„ 

o 

g, 

o 

g, 

o 


o, 

o 

g, 

o 


g, 

o 


end 


g, 

o 

g, 

o 

g, 

o 

g, 

o 


end 


o,  o, 
o  o 


end 


disp ( ' Temp : ' ) , disp (temp) 
disp('H  Vel :'), disp (velocity ) 
disp('V  Vel disp (vVel ) 
di sp ( ' Row :  ' )  ,  di sp ( r ) 
disp ( ' Col : ' ) , disp ( c ) 

disp ( ' Strain  Rate ' ) , disp ( srData (r , c ) ) 
error (' Strain  Rate  Out  of  Range') 


%  EVALUATE  STRAIN  AT  MAX  STRESS  FAILURE  AREA 


o. 

o 

o, 

o 

g, 

o 

g, 

o 

g, 

o 

o, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 

g, 

o 


failureSMS  =  zeros (CycleNum, TracerPoints ) ; 

f ailureSumSMS  =  zeros (CycleNum, 1 ) ; 

A  =  2 . 24700e-2 ; 

B  =  -5.5160e-2; 

C  =  6 . 04400e-3 ; 

f ailureCritSMS  =  A* ( srData . AB)  +  C; 

for  row=l : r 

for  col=l : c 

if  row>l  &&  failureSMS (row-1 , col ) ==1 
failureSMS (row,  col ) =1 ; 

end 

if  strainData (row,  col) >=f ailureCritSMS (row, col) 
failureSMS (row, col ) =1 ; 

end 

end 

end 

failureSMS  =  failureSMS . *vf Data;  %unitless 

for  iter  =  1 : length ( f ailureSumSMS ) 

f ailureSumSMS ( iter , 1 )  =  sum ( failureSMS ( iter ,:)) ;  %unitless 

end 

damAreaSMS  =  f ailureSumSMS*meshSize;  %mmA2 


for  i  =  2: CycleNum 

distanceSlid  =  velocity*time ( i ) *1000 ; 

%WR_SMS ( i )  =  100*damAreaSMS (i) /distanceSlid; 

end 


SMS_Average  =  100*damAreaSMS (end) /distanceSlid 
WR_SMS  =  100*damAreaSMS/distanceSlid;  %mmA2/mm 
disp ('Strain  at  Max  Stress  Failure  Mechanism  Evaluated...') 

%%  Evaluate  JFrac  Data  -  method  2 
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final_cells  =  JFData ( [end]  ,  : ) ; 
failureJF  =  final_cells . *vfData (101,  :  )  ; 
failureSumJF  =  sum ( failureJF ) ; 
damAreaJF  =  f ailureSumJF*meshSize ; 

WR_JF  =  100*damAreaJF (end) /distanceSlid; 

%for  i  =  2:CycleNum 

%  distanceSlid  =  velocity*time (i) *1000; 

%  WR_JF(i)  =  100*damAreaSMS (i) /distanceSlid; 

%end 

JF_Average  =  100*damAreaJF (end) /distanceSlid 
disp ( ' Johnson-Cook  Fracture  Criterion  Evaluated.. 

9-  9- 
o  o 


%Calculate  Wear  at  each  interval 
plot (time,  WR_JF) 
xlabel('Time  (s)') 

ylabel('Wear  Rate  per  Unit  Width  (mmA2/mm) ' ) 
title ( ' 2um  Asperity  Wear') 

figure 

plot (time, damAreaJF) 
xlabel ( ' Time  ( s ) ' ) 

ylabel('Area  of  Damaged  Elements (mmA2 )' ) 
title ( ' 2um  Asperity  Damage  Area') 
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Appendix  C.  Wear  Damage  Area  vs.  Time  Plots 


Below  are  all  of  the  figures  generated  from  looking  up  the  maximum  wear 
damage  over  the  simulation  time.  At  the  end  of  the  section  the  maximum  values  are 
included  tabularly,  along  with  the  sliding  distance  that  each  maximum  occurred  at. 
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Wear  per  unit  width  (mmA2/mm)* 

Strain  at  Max  Stress 

vonMises 

Johnson  and  Cook  Fracture 

200mps 

500mps 

lOOOmps 

1500mps 

200mps 

500mps 

lOOOmps 

1500mps 

200mps 

500mps 

lOOOmps 

1500mps 

6micron 

0.000243 

0.000146 

0.000383 

0.000599 

0.000124 

0.000134 

0.000299 

0.000395 

5.85E-05 

5.21E-05 

0.000176 

0.000278 

4micron 

6.26E-05 

8.21E-05 

1.77E-05 

0.00031 

2.09E-05 

5.27E-05 

1.83E-05 

0.000294 

5.55E-05 

4.73E-05 

7.71E-05 

0.000126 

2micron 

3.22E-05 

2.76E-05 

6.61E-05 

5.41E-05 

7.46E-06 

1.44E-05 

4.1E-05 

7.11E-05 

7.98E-06 

7.83E-06 

1.89E-05 

1.29E-05 

lmicron 

2.32E-05 

9.91E-06 

2.71E-05 

1.42E-05 

1.93E-06 

1.57E-06 

6.98E-06 

1.41E-05 

4.67E-07 

7.47E-07 

2.67E-06 

5.28E-06 

0.5micro 

n 

8.82E-06 

8.65E-06 

3.88E-06 

5.6E-06 

2.44E-06 

7.43E-07 

3.67E-07 

9.18E-07 

0 

2.36E-08 

1.98E-08 

3.56E-07 

Distance  slid  at  maximum  wear  (mm 

* 

Strain  at  Max  Stress 

vonMises 

Johnson  and  Cook  Fracture 

200mps 

500mps 

lOOOmps 

1500mps 

200mps 

500mps 

lOOOmps 

1500mps 

200mps 

500mps 

lOOOmps 

1500mps 

6micron 

0.0216 

0.01965 

0.0364 

0.04965 

0.02 

0.02015 

0.00891 

0.00996 

0.0198 

0.014 

0.0364 

0.04965 

4micron 

0.00872 

0.0183 

0.00202 

0.01845 

0.00872 

0.0033 

0.00202 

0.000293 

0.01636 

0.022 

0.0247 

0.0243 

2micron 

0.01186 

0.01015 

0.0218 

0.0051 

0.0124 

0.00715 

0.0153 

0.005025 

0.0088 

0.01015 

0.0201 

0.0051 

lmicron 

0.011 

0.011 

0.0107 

0.002895 

0.011 

0.011 

2.29E-03 

0.002895 

0.00958 

0.0079 

1.31E+04 

0.010605 

0.5micron 

0.0055 

0.0055 

0.0055 

0.005505 

0.0055 

0.0055 

0.0055 

0.005505 

0.00E+00 

0.00055 

0.00055 

0.00513 
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Appendix  D.  Single  Asperity  Wear  Rates  and  Wear  Volume  Calculations 


The  following  tables  contain  the  wear  rates  generated  from  CTH  that  describe  the 
maximum  wear  rates  at  200,  500,  1000,  and  1500  meters  per  second  for  asperities  with  6, 
4,  2,  1,  and  0.5  radii.  These  wear  rates  are  then  used  to  calculate  single  asperity  wear  rates 
and  averaged  together  to  calculate  spherical  correction  factors.  These  calculations  are 
also  shown  in  the  tables  in  this  section. 


Table  21:  06micron  Wear  Rates  using  Johnson  and  Cook  Failure  Criterion 


Wear  rate  per  unit  width  (Wuw) 

200 

500 

1000 

1500 

Distance 

06 

microns 

0.002952 

0.003721 

0.004844 

0.005605 

0.00 

sqrt(6A2- 

6A2) 

04 

microns 

0.003395 

0.002148 

0.003122 

0.005203 

4.47 

i 

CM 

< 

CD^ 

tr  cm 

O"  < 

CO 

02 

microns 

0.000907 

0.000771 

0.00094 

0.002536 

5.66 

sqrt(6A2- 

2A2) 

0 

microns 

0 

0 

0 

0 

6 

Table  22:  Area  under  the  Curve  Calculations  for  6micron  Wear  Rates  using  Johnson  Cook  Failure 

Criterion 


Area 

200mps 

500mps 

lOOOmps 

1500mps 

0.000988709 

0.003515341 

0.003847749 

0.000898107 

1.52E-02 

9.60E-03 

1.40E-02 

2.33E-02 

0.001479842 

0.000819275 

0.00129852 

0.001586831 

1.08E-03 

9.18E-04 

1.12E-03 

3.02E-03 

0.000154255 

0.000131099 

0.000159772 

0.000431165 

Total 

area 

0.016898643 

0.014985418 

0.02038077 

0.029192578 

x2 

0.033797286 

0.029970836 

0.04076154 

0.058385155 
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Table  23:  Spherical  Correction  Factor  Calculations  for  6micron  using  Johnson  and  Cook  Failure 

Criterion 


m/s 

Single 
asperity 
wear  rate 
(Wsa) 
mmA3/mm 

Ratio 

(Wsa/Wuw) 

um 

200 

3.380E-02 

1.14E+01 

500 

2.997E-02 

8.05E+00 

1000 

4.076E-02 

8.42E+00 

1500 

5.839E-02 

1.04E+01 

9.58E+00 

Table  24:  6micron  Wear  Rates  using  the  Strain  at  Max  Stress  Criterion 


Wear  rate  per  unit  width  (Wuw) 

200 

500 

1000 

1500 

Distance 

06 

microns 

0.011247 

0.007437 

0.010513 

0.012055 

0.00 

sqrt(6A2- 

6A2) 

04 

microns 

0.007182 

0.004489 

0.008758 

0.016797 

4.47 

i 

CM 

< 

CD^ 

tr  CM 
O"  < 

CO 

02 

microns 

0.002711 

0.002715 

0.003031 

0.010604 

5.66 

sqrt(6A2- 

2A2) 

0 

microns 

0 

0 

0 

0 

6 

Table  25:  Area  under  the  Curve  Calculations  for  6micron  Wear  Rates  using  Strain  at  Max  Stress 

Criterion 


200mps 

500mps 

lOOOmps 

1500mps 

0.009085256 

0.006588593 

0.003922881 

0.010598313 

3.21  E-02 

2.01  E-02 

3.91  E-02 

7.51  E-02 

0.002659861 

0.001055453 

0.003407672 

0.00368475 

3.23E-03 

3.23E-03 

3.61  E-03 

1.26E-02 

0.000460915 

0.000461502 

0.000515226 

0.001802663 

Total 

area 

0.047534267 

0.03140004 

0.05060028 

0.08258925 

x2 

0.095068535 

0.06280008 

0.101200561 

0.1651785 
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Table  26:  Spherical  Correction  Factor  Calculations  for  6micron  using  Strain  at  Max  Stress  Failure 

Criterion 


m/s 

Single 
asperity 
wear  rate 
(Wsa) 
mmA3/mm 

Ratio 

(Wsa/Wuw) 

um 

200 

9.507E-02 

8.45E+00 

500 

6.280E-02 

8.44E+00 

1000 

1.012E-01 

9.63E+00 

1500 

1.652E-01 

1.37E+01 

1.01E+01 

Table  27:  2micron  Wear  Rates  using  Johnson  and  Cook  Failure  Criterion 


Wear  rate  per  unit  width  (Wuw) 

200 

500 

1000 

1500 

Distance 

2 

microns 

0.000907 

0.000771 

0.000940 

0.002536 

0.00 

sqrt(2A2- 

2A2) 

1  micron 

0.000049 

0.000095 

0.000000 

0.000498 

1.73 

sqrt(2A2- 
1 A2) 

0.5 

microns 

0.000000 

0.000043 

0.000036 

0.000069 

1.93 

sqrt(2A2- 

0.5A2) 

0 

microns 

0 

0 

0 

0 

2 

sqrt(2A2-0) 

Table  28:  Area  under  the  Curve  Calculations  for  2micron  Wear  Rates  using  Johnson  Cook  Failure 

Criterion 


Area 

200mps 

500mps 

lOOOmps 

1500mps 

0.000743604 

0.00058599 

0.000813896 

0.001765468 

8.44E-05 

1.64E-04 

3.54E-1 0 

8.62E-04 

4.82306E-06 

5.1 0522E-06 

-3.56493E- 

06 

4.2401 9E-05 

0.00E+00 

8.50E-06 

7.13E-06 

1.37E-05 

0 

1.50285E-06 

1 .26034E-06 

2.42596E-06 

Total 

area 

0.000832806 

0.000764784 

0.000818722 

0.002685889 

x2 

0.001665613 

0.001529569 

0.001637444 

0.005371778 
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Table  29:  Spherical  Correction  Factor  Calculations  for  2micron  using  Johnson  and  Cook  Failure 

Criterion 


m/s 

Single 
asperity 
wear  rate 
(Wsa) 
mmA3/mm 

Ratio 

(Wsa/Wuw) 

um 

200 

1.666E-03 

1.84E+00 

500 

1 .530E-03 

1.98E+00 

1000 

1 .637E-03 

1.74E+00 

1500 

5.372E-03 

2.12E+00 

1.92E+00 

Table  30:2micron  Wear  Rates  using  the  Strain  at  Max  Stress  Criterion 


Wear  rate  per  unit  width  (Wuw) 

200 

500 

1000 

1500 

Distance 

2 

microns 

0.002711 

0.002715 

0.003031 

0.010604 

0.00 

sqrt(2A2- 

2A2) 

1  micron 

0.002113 

0.000901 

0.002536 

0.004896 

1.73 

sqrt(2A2- 
1 A2) 

0.5 

microns 

0.001604 

0.001573 

0.000706 

0.001017 

1.93 

sqrt(2A2- 

0.5A2) 

Omicrons 

0 

0 

0 

0 

2 

sqrt(2A2-0) 

Table  31:  Area  under  the  Curve  Calculations  for  2micron  Wear  Rates  using  Strain  at  Max  Stress 

Failure  Criterion 


Area 

200mps 

500mps 

lOOOmps 

1500mps 

0.000518283 

0.001570781 

0.000428277 

0.004942856 

3.66E-03 

1.56E-03 

4.39E-03 

8.48E-03 

5.04025E-05 

-6.65686E- 

05 

0.000181203 

0.000384037 

3.18E-04 

3.12E-04 

1.40E-04 

2.01  E-04 

5.61284E-05 

5.50652E-05 

2.47053E-05 

3.55967E-05 

Total 

area 

0.004601684 

0.003431116 

0.00516664 

0.014044108 

x2 

0.009203367 

0.006862233 

0.01033328 

0.028088217 
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Table  32:  Spherical  Correction  Factor  Calculations  for  2micron  using  Strain  at  Max  Stress  Failure 

Criterion 


m/s 

Single 
asperity 
wear  rate 
(Wsa) 
mmA3/mm 

Ratio 

(Wsa/Wuw) 

um 

200 

9.203E-03 

3.39E+00 

500 

6.862E-03 

2.53E+00 

1000 

1 .033E-02 

3.41  E+00 

1500 

2.809E-02 

2.65E+00 

3.00E+00 
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